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ABSTRACT 

Cosmological parameter measurements from CMB experiments such as Planck, ACTpol, SPTpol 
and other high resolution follow-ons fundamentally rely on the accuracy of the assumed recombination 
model, or one with well prescribed uncertainties. Deviations from the standard recombination history 
might suggest new particle physics or modified atomic physics. Here we treat possible perturbative 
fluctuations in the free electron fraction, Xc{z), by a semi-blind expansion in densely-packed modes 
in redshift. From these we construct parameter eigenmodes, which we rank order so that the lowest 
modes provide the most power to probe the Xc{z) with CMB measurements. Since the eigenmodes are 
effectively weighed by the fiducial history, they are localized around the differential visibility peak, 
allowing for an excellent probe of hydrogen recombination, but a weaker probe of the higher redshift 
helium recombination and the lower redshift highly neutral freeze-out tail. We use an information- 
based criterion to truncate the mode hierarchy, and show that with even a few modes the method 
goes a long way towards morphing a fiducial older Recfast Xej(z) into the new and improved 
CosmoRec and HyRec X^jiz) in the hydrogen recombination regime, though not well in the helium 
regime. Without such a correction, the derived cosmic parameters are biased. We discuss an iterative 
approach for updating the eigenmodes to further hone in on Xc,f(z) if large deviations are indeed 
found. We also introduce control parameters that downweight the attention on the visibility peak 
structure, e.g., focusing the eigenmode probes more strongly on the Xe(z) freeze-out tail, as would be 
appropriate when looking for the X^ signature of annihilating or decaying elementary particles. 



1. INTRODUCTION 

The Planck Surveyoi[^is now well into its mission, ob- 
serving the temperature and polarization anisotropics of 
the cosmic micro wave background (CMB) with unprece- 
dented accuracy (Planck HFI Core Team et al. 2011[ 
mm. Both AC'l' (e.g. see|Hajlan et al! 

Das et all [2011] ) and SP'l' 



Mennella et al.ll20l^ 

^unklcy ct a 



2010 



( [Lueker et al.i2010f Vanderlihde et al.|2Ulu| ) are pushing 



the frontier or CmB power s pectra at small scales 
and in the near f uture SPTPOlP] jMcMahon et al.|2009 1 
and ACTPOI0 ( [Niemack et~al||2010[ ) will provide ad- 
ditional small scale A'-mode polarization data, comple- 
menting the polarization power spectra obtained with 
Planck and further increasing the significance of TT 
power spectra. 

Using these datasets, cosmologists will be able to de- 
termine tlie_^e;^_cosm£logicaJ__gaxM with high pre- 
cision (The Planck Collaboration"2006', Tauber et al. 
2010 1 , makmg it possible to distinguis h between v ari- 
ous models of inflation (e.g. see Komatsu et al.|2011 for 
recent constraints from WMAP) by measuring the pre- 
cise value of the spectral index of scalar perturbations, 
rig, and constraining its possible running, rirun, as well 
as the tensor-to-scalar ratio, r. In addition many non- 
standard extensions of the minimal inflationary model 
are under discussion, a nd the observability of these pos- 
sibilities with PLANCK(|The Planck Collaboration||2006 1 
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and future CMB experiment is being considered. 

These encouraging observational prospects have moti 
vated vario us inde pendent g roups (e.g. see 'Dubrovich 



fc Grachev 2005; C hluba fc S unyacv 2006b; Kholupcnko 
& Ivanchik 2006; S witzer fc Hhata 2008; Wong fc Scott 



2007; Rubiho-Marti'n et al.|2008 KarshenEoiin & Ivanovl 
2008; Hirata 2008; Chluba £~5unyaev 2008rJcntschu^ 
2009. .Labzowsky et al.,.2009| |Grin "S: Hirata„2010, .Ah , 
Haimoud & Hirata||2010p to assess how uncertainties in 



the theoretical treatment of the cosmological recombina- 
tion process could affect the science return of Planck 
and future CMB experiments. The precise evolution of 
the free electron fraction, X^, with time influences the 
shape and position of the peak of the Thomson visi bility 
function, which defines the last scattering su rface ( Sun- 



^ http : //www .physics .princeton. edu/ act/ 



yaev fc Zeldovich|[l970| [Peebles fc Yu||1970| , and hence 
controls how photons and baryons decouple as electrons 
recombine to form neutral helium and hydrogen atoms. 
Consequently, the ionization history changes the acoustic 
oscillations in the photon-baryon fluid during recombina- 
tion and therefore directly affects the CMB temperature 
and polarization power spectra. For the analysis of fu- 
ture CMB data this implies that in particular close to 
z ^ 1100 the ionization history better be understood at 
the - 0.1% level. 

Probing the ionization history in time is equivalent to 
probing it in space with the light cone relating the two. 
Thus what we try to do in this paper, namely to come up 
with optimized probing functions for the recombination 
history, is quite akin to creating probes of the spatial 
structure of the boundary between HII and neutral hy- 
drogen regions. Here of course we look from neutral to 
ionized, the cosmological recombination problem being 
an inside-out HII region, except in a predominantly elec- 
tron scattering regime with a very large photon to baryon 
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ratio which lowers the transition temperature between 
ionized and neutral. 

The old rec ombination standard was set by Recfast 
( [Seager et a"L] [l999, 2000), but its reUability for the pre- 



cision cosniology was brought into question, e.g., by |Sel- 
jak et al. ( |2003[ ) . For the standard six parameter cos- 



moiogy in particular our ability to measure the precise 
value of Us and the baryon content of our Universe may 
be compromised if modifications to the recombination 



model of Recf ast are neglected ( Rubiho-Martm et al. 
2010 Shaw fc Chluba 2011) , introducing biases ot a tew 
a tor FLANCK. ' ' 

Currently it appears that all important corrections 
to the standard recombinatio n scenario (SRS hereafter) 



have been identified (e.g., see Fendt et al.||2009[ [Rubiiio- 
Martm et al.||20l"0 for an overview). I'he new reco mbi- 
nation co des, UOSMoRec ( Chluba fc Thomas|2011 ") and 
HyRec ( Ali-Hai'moud fc Hirata, 2011] ) both account for 
these modihcations to the SRS, superseding the physical 
model of Recfast and allowing fast and accurate com- 
putation of the ionization history on a model-by-modcl 
basis. CosmoRec and HyRec presently agree at a level 
of ~ 0.1% — 0.2% during hydrogen recombination, so that 
from standard recombination physics little room for big 
surprises seems to be left. 

However, what if something non-standard happened? 
What if something was overlooked in the standard recom- 
bination scenario? From the scientific point of view the 
ionization history is a theoretical ingredient to the cos- 
mological model, which usually is assumed to be precisely 
known and not subject to direct measurement. Clearly, 
it is important to estimate the possible level of uncer- 
tainty in the recombination model and to confront our 
understanding of the recombination problem with direct 
observational evidence. Here we describe how well future 
cosmological data alone are able to constrain possible de- 
viations from the SRS. 

In the past, several non-standard extensions of the re- 
combination scenario have been considered. These in- 
clude models of delayed recombination, in which hypo- 
thetical sources of extra photons that can lead to ioniza- 
tions or excitations of atoms are introduc ed using sim- 
ple parametrizations (Peebles e t al. 2000). In particu- 



lar, models of deca,yinQ [e.^., seelchen fe 
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Zhang et al.|[2007 ) and ann ihilatin g particles {e.g.' 
Fadma nabhan fc Fi nkbeinerf 2005i Zh ang et al.|2006" 
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fundamental constants (e.g., see Kaplinghat et al. 


1999 


Scoccola et al. 


2009 


Galli et al. 


2009f could atlect the 



All these ideas rely on a specific model for the (phys- 
ical) process under consideration, with the derived con- 
straints depending on the chosen parametrization. This 
minimizes the number of additional parameters, but does 
not allow us to answer questions about more general per- 
turbations around the SRS and how well they can actu- 
ally be constrained. 

Here we approach this problem in a different way. We 
introduce perturbations to the SRS over a wide range 
of redshifts around hydrogen (z ~ 1100) and helium 
(z ~ 1800) recombination, using different basis functions. 
We then compute the corresponding signals in the CMB 
power spectra and perform a principle component de- 



composition to obtain eigenmode functions, ordered with 
respect to the level at which they can be constrained by 
the data. We study in detail how the eigenmodes depend 
on different experimental settings, the fiducial model, as 
well as the chosen parametrization for the recombination 
perturbations. 

Our metho d is similar to the one used by |Mortonson| 
fc Hu (20081, where the eigenmodes for different reion- 
ization scenarios (6 < z < 30) were constructed. How- 
ever, here we explicitly construct the mode functions at 
redshifts z > 200, with particular attention to the de- 
pendence of the eigenmodes on different assumptions. 
We investigate how to use our prior knowledge of possi- 
ble perturbations of the ionization history to choose the 
parametrization which is more preferred by the data. We 
also carry out a careful convergence study and show the 
equivalence of different basis functions (e.g., triangles, 
Gaussian bumps, Fourier series and Chebyshev polyno- 
mials). We particularly focus on the helium recombina- 
tion problem, showing that in the absence of very tight 
constraints on the hydrogen recombination, we will not 
be able to unravel well remaining uncertainties in helium 
recombination with CMB data. 

Similarly, small changes in the freeze-out tail of recom- 
bination are only weakly constrained, if possible ambigu- 
ities during hydrogen recombination are included. 

Details of the general methodology to construct the 
eigenmodes for perturbations to ionization history are 
given in § [2] In § |3] we computed different eigenmodes 
over a rather wide redshift range (z e [200,3000]). At 
the end of § [3j we develop a criterion which allows us to 
truncate the Hierarchy of the eigenmodes based on their 
information content. In § |4] the modes are applied to 
two specific examples of ionization scenarios, illustrating 
how the method should be used with real CMB data. 
At the end of this section, we also discuss how the ap- 
proach should be iterated if hints toward a considerable 
difference between the assumed and true model of recom- 
bination are indicated by the data. We close the paper 
by a brief discussion. 

2. METHODOLOGY 

In this section we introduce the approach and 
parametrization used to construct the principle compo- 
nents, or the eigenmodes, which will be used to describe 
possible corrections to the recombination scenario. Our 
method is mainly driven by the assumption of small rel- 
ative perturbations arou nd the fiducial model co mputed 
with the Re cfast code ( [Seager et al.| ( |1999| ; see [Wong 
et al. ( 2008 ) for recent updates). As an example we have 
m mmd the recombination co rrections obtained with re- 
fined recombination models (Chluba fc Thomas] |201l| 
Ali-Haimoud fc Hirata 2011). However we also briefiy 
discuss the possibility to constrain significant changes 
in the freeze-out tail of recombination and modes that 
mainly focus on helium recombination. 

Throughout this paper the cosmic parameters, which 
will be referred to as the standard (cosmological) param- 
eters, are {flhh^ ,i^dmh^ , Ho,T,ns, Ag) as measured by 
WMAPT]^ unless stated otherwise. In several cases we 
also vary Yp as a seventh parameter. Lensing is included 

^ http:/ /lambda. gsfc.nasa.gov/product/map/dr4/params 
/lcdm_sz_lens_wmap7.cfm 
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in all simulations if not explicitly stated otherwise. 



2.1. The standard recombination scenario 

As explained in the introduction, the cosmological re- 
combination history is one of the major theoretical in- 
puts for computations of the CMB anisotropics. Conse- 
quently, high precision unbiased cosmic parameter mea- 
surements from current and future CMB experiments re- 
quire a sufficiently accurate model for hydrogen and he- 
lium recombination. 

The ionization fraction for the SRS is shown in the left 
panel of Fig. [l] It was calculated using Recfast vl.4.2, 
which accounts for some of the modification to helium 



recom bination (Kholup enko et a l. 2007; Swi tzer fc Hi- 
rata|2008j Rubino-Martin et al.l2008,,Chluba & Sunyaev 

bij 



2010p using fudge parameters, but neglects detailed 7a^ 



diative transfer corrections (see [Chluba fc Thomas|2011 



Ali-Haimoud 
around z 



fc Hirata 
TTUU 



2011 



and references therein) 
curve corresponds to an 



ionization fraction with the meas ured temperature of the 
CMB radiation, Tcmb 2.726K ([Fbcscn 2009). For com- 
parison and to illustrate the temperature dependence of 
the ionization history, the ionization fraction correspond- 
ing to Tcmb = 3K is also plotted (dashed curve). A 
larger value of Tcmb means more photons in the Wien 
tail of the CMB blackbody, so that the matter is kept 
ionized until lower redshift. 

On the right the corresponding differential visibility 
functions (or visibility functions for short) are plotted: 



dry 



(1) 



where rj is the conformal time and r is the Thomson 
scattering optical depth from redshift z to now. 

The visibility function describes the probability that a 
photon we observe today last scattered off free electrons 
at a certain position along the line of sight. The CMB 
anisotropics formed mainly during the epoch of hydrogen 
recombination defined by the peak of the visibility func- 
tion located at redshift z ~ 1100. They are thus most 
sensitive to changes around the maximum of visibility. 
For example, an increase in the width of the visibility 
bump corresponds to a more extended or slower recom- 
bination process, leading to more Thomson scatterings 
of photons off free electrons. These scatterings lead to 
the cancellation of the CMB anisotropics along the line of 
sight on scales comparable and smaller than the recombi- 
nation width, while enhancing the polarization signal on 
larger scales. The location of the maximum of the visibil- 
ity function for an assumed cosmological model, on the 
other hand, determines the distance to the last scatter- 
ing surface. This in turn determines the positions of the 
peaks of the CMB power spectra. Similarly, any change 
in the ionization history, through affecting the visibility, 
would lead to (possibly measurable) changes in the CMB 
power spectra. 

As the right panel of Fig. [T] indicates, at high redshifts 
z > 1400 the visibility function falls off very quickly. At 
those times the number of free electrons is still so large 
that scatterings occur very frequently and the mean free 
path is very short. Consequently, the part of the ioniza- 
tion history which is connected to helium recombination 



mainly affects the damping tail of the CMB anisotropics, 
but even there the effect is rather moderate, as in the red- 
shift range 

n 1400 < z< 3000 helium can at most alter 
the number of electrons by ~ 8%. 

2.2. Choice of perturbation parametrization 

We now introduce small perturbations to the SRS that 
allow construction of the eigenmodes. 

There are different ways to parametrize the pertur- 
bations to the ionization history in a (semi-)model- 
independent way. For example, to study how well the 
low redshift ionization history {z in the i nterval [6, 30]) 
can be constr ained by future CMB data, |Hu fc Holder] 
( 2003 1 and ,Mortonson fc Hu ( 2008 1 used the changes in 
the ionization traction in ditterent redshift bins, dXc{z) — 
const, to parametrize the uncertainties. This is a valid 
choice for the low redshift region, because our ignorance 
of the underlying model of reionization does not sug- 
gest any preferred non-uniform weighting of the per- 
turbations at different redshifts. In this regime dXc{z) 
probes the ionization fraction itself and not perturba- 
tions guided by a fiducial model. The results from this 
choice of parametrization are shown to be fiducial model- 
independent which is expected due to the weak signal 
from the reionization process. 

In contrast to this, at high redshifts (z ~ 1100) there 
is strong theoretical support for the exhaustively studied 
model of recombination in the realm of standard atomic 
physics and radiative processes. Also as it was mentioned 
in the introduction, the current generation of CMB data 
is sensitive to changes in at the level of a few percent. 
Therefore the main assumption in this paper is that the 
fiducial model for the ionization history, X^'^{z), is close 
to the true underlying history, Xe{z), which we are look- 
ing for. We call this method semi-blind emphasizing our 
belief in the SRS as the framework of recombination, with 
the search for deviations being limited to small perturba- 
tions around this reference model. The goal is to detect 
or place upper limits on these potential small deviations. 
Clearly, if data point toward significant deviations from 
the SRS, an itera tive approach should be adopted, as will 
be discussed in § |4.4[ 

With small deviations in mind we can write: 

X,iz)^X^^iz) + SX,{z), 

with jfJATcl/X^'^ <C 1. A natural parameter to describe 
the perturbation is then the relative deviation in the ion- 
ization fraction: 

Su{z)=dXc{z)/X^'^{z) with |(5m(2)| < 1. (2) 

This parametrization has the advantage of always satis- 
fying the necessary condition X^ > 0. It is also straight- 
forward to fulfill the Xc < A"o,max condition in the simu- 
lations, where Xe,max is determined by Yp, the primordial 
helium mass abundance, through X^ max — 1 -|- yp/2(l — 
Yp). The parametrization in Eq. ([2]) weights possible 
perturbations at different redshifts oy the fiducial ion- 
ization fraction. This implies that for du{z) = const 
the absolute change in the ionization fraction \6Xe\ is 
down- weighted in the freeze-out tail of Xc (z < 800; see 

^ The recombination of doubly ionized helium ends around red- 
shift z ~ 5000. 



4 



Far hang et al. 




2000 



Fig. 1. — The cosmological ionization histo ry (left). Xp = Nc/[Np 
recombination scenario with TcMB = 2.726K | |Fixsen|20 0"9| contrasted 
density of ionized and neutral hydrogen, while JVe denotes the number 

Fig. T]) , compared to perturbations around maximum vis- 
ibility {z ~ 1100) where X, ~ 1. Throughout this paper 
Su{z) as defined in Eq. ([2| will be our main choice of 
parametrization. 

A more general parameter which includes the above 
parametrization as a special case is given by: 



Su{z) = dX,{z)/[X^''iz) + a{z)] 



(3) 



where 17(2:) > can be a constant or otherwise convenient 
function of redshift allowing to focus on different redshift 
ranges of interest. In particular, when considering possi- 
ble modifications to the ionization history introduced by 
energy injection, e.g. because of annihilating dark mat- 
ter, or dec aying r elic partic les (Chen & Kamion kowski 
Padman abhan fc Finkbeincr, ,2005,; Zhang et al. 
2007; H fitsi et al.||2009| [Slatyer et ai.||2009| [Htitsi 



2004 



2006 



et al.|2011| ) , where the freeze-out tail of recombination is 
disturbed, a value of a X^'^ would be a proper choice, 
giving higher weigh t to the pert urb ations in the lower 
redshift part (see § 3.4 and Fig. 10). In the limit of a 
high value of a relative to the fiducial X^ the parame- 
ters approach Su{z) = SXc{z) which uniformly weights 
perturbations at different redshifts. This, as already dis- 
cussed, is a good choice for regions where there is no 
strong a priori belief in the underlying model or if the 
redshifts of interest have comparatively low X^ where 
6u{z) with a — does not lead to strong enough signals 
to probe. In principle, a conveniently chosen redshift de- 
pendent cr(z) is a tool to effectively incorporate our prior 
knowledge of the ionization history in the parametriza- 
tion of its perturbations. For example, with Su{z) de- 
fined by Eq. (|3| one can smoothly interpolate between 
relative and absolute perturbations to X^, at high and 
low redshifts respectively. Also it is clear that one can 
focus on different parts of the recombination history by 
limiting the redshift range over which the eigenmodes 
are constructed [e.g., just on reionization (0 < 2: < 30) 
or helium recombination (1400 < 2 < 3000)]. 

2.2.1. Alternative parametrizations 



+ -/Vjji], and differential visibility function (right) for the standard 
to a case with Tcmb = 3K. Here A'^p and A^hi represent the number 
density of free electrons. 

Instead of directly perturbing the ionization fraction, 
as we chose here, it is plausible to parametrize possi- 
ble changes in the physical sources of perturbation to 
the ionization history, such as energy injection in the 
medium which leads to excitation or ionization of atoms, 
or the Lya escape probabili ty during recombin ation (see 
introduction). For example Mitra et al. (2010) chose the 
number of photons in the iUM per baryon m collapsed 
objects as the parameter to study the low redshift ioniza- 
tion history. Alternatively, one could modify the fudge 
factors or functions in Recfast, or alter the expansion 
rate given by the Hubble factor, H{z). 

Each of these possibilities implies different priors on 
the regions that can be altered and, e.g., in the case of 
H{z) other aspects of the cosmological model are also 
affected. They also cover, in general, only a limited class 
of changes to the recombination history. When inter- 
ested in perturbations to the ionization history, Xc is 
the physical quantity which, via the visibility function 
and the optical depth, most directly enters the Boltz- 
mann equations which describe radiation anisotropics 
and ca n be solved using the Boltzma nn codes such as 
Camb ( Lewis et al.|200() ) or Cmbfast ( [Seljak fc Zaldar 



riaga||1996 ]. I'he ionization fraction has the additional 
advantage, over the visibility and the optical depth, of 
being straightforward to limit to physically allowed val- 
ues. The nearly direct mathematical encounter of X^ 
with CMB anisotropics guarantees that any perturbation 
in the plasma that would lead to changes in the radiation 
anisotropics should go through and thus be reflected in 
Xc- Therefore the relative changes in Xc constitute our 
preferred physical parameters. 

We close by mentioning that, it is also theoretically 
possible to consider different variables for time such as 
(conformal) time, optical depth and scale factor. How- 
ever, for our purpose we choose to work with redshift 
to describe temporal dependence. In principle different 
parametrizations, if they cover the same range of phys- 
ical perturbations, can be transformed to one another 
with the proper change of the a priori distribution of pa- 
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rameters. Here, in the absence of physicaUy motivated 
constraints, a uniform prior is assumed for perturbations 
at different redshifts regardless of parametrization (here, 
e.g., for various values of a in Eq. ([s])). If the pertur- 
bation is strongly constrained by data, the choice of the 
prior would not play a major role. 

2.3. Basis Junctions and their different characteristics 

Having chosen the parametrization, we now expand 
the perturbations in a discrete set of mode (or basis) 
functions, ifi(z): 



N 



^min ^ ^ ^ ^inax (4) 



and 6u{z) — elsewhere. Here r[z) is the residual and 
yfs are the parameters defining the strength on the mode 
ipi{z). Often we take ipi{z) to be localized in z about a 
knot value Zi , but this is not necessary. We can, for ex- 
ample, choose the ^Piiz) to form a complete orthonormal 
set in which case N ^ oo and the residual r would be 
zero. Below, we will discuss different possibilities for the 
choice of the mode functions. 

2.3.1. Localized basis functions 

We first investigated localized Gaussian and triangular 
bumps as mode functions. Both can be considered as ap- 
proximations to the Dirac (5-function. We define the ith 
basis function centered at redshift Zi and having width 
CTj by: ^ 

tpi{z) cx exp (-^—^^^j (5) 



for the Gaussian case and by 



ip,{z) cx 







\z - Zj| < (Ji 

otherwise, 



(6) 



for the triangles. Triangular bumps were used earlier 
in the principle c omponent analysis of different reion- 
ization scenarios (|Hu & Holder 2003; Mortonson & Hu 
20081. In some circumstances, tlie sharp edges in the 



triangles could cause numerical problems. Smoothed lo- 
calized functions such as Gaussians and the M4 splines 
introduced below have numerical advantage. 

As basis functions, it is more convenient if the set of 
<y9i's is an orthogonal set. For this, there should be no 
overlap between different bumps. On the other hand, 
there is no way to cover the whole redshift range -a 
necessary condition for completeness- with a finite num- 
ber of non-overlapping bumps. However, depending on 
the problem of interest, the width and separation of the 
(overlapping) bumps can be properly chosen to ensure all 
points in the redshift interval have been covered, while at 
the same time the orthogonality is not strongly violated. 

Instead of Gaussian and triangular bumps, one can also 
adopt an approach similar to that used in Smoothed- 
Particle Hydrodynamics (SPH), and think of the basis 
functions as window functions (or kernels) used to in- 
terpolate the properties of particles to any point in the 
medium. For us, the par ticles would be the spline knots 
(e.g. see |De Boor||2001|) at the specific with the as- 



sociated magnitude yi. There is a smoothing length h 
associated with the kernel over which the properties of 



the particles are smoothed. A commonly used kernel 
(other than th e Gaussian consid ered above) is the cubic 
Mi spline (e.g. Monaghaii||2005 ) , defined by: 



ipi{z) cx Mi{\z - z,\) 



0<-z< 1; 
1 <g<2; 
q>2- 



(7) 



where q = \z — Zi\/h. Whereas the Gaussian kernel has 
non-zero contributions from every redshift (though the 
range is usually truncated beyond about 3(t), the cubic 
spline is compact , rea ching zero for particles beyond 2h. 
(As discussed in § 3.2 we have used the smoothing length 
h — l.bSz where dz is the particle separation). 

We modified the publicly available code CAM^to sim- 
ulate CMB power spectra for a more general recombina- 
tion scenario that includes perturbations on top of the 
SRS. Introducing narrow features into the ionization his- 
tory also required an increase in the redshift sampling of 
A"o. We checked the numerical convergence and stability 
of the results by using high accuracy settings. 

As examples for localized perturbations, Fig. [2] shows 
three perturbation functions 6u{z) = 51n(Xe) based on 
M4 splines (left panel) and the corresponding Ci response 
in TT and EE. The perturbations are located at differ- 
ent redshifts and have equal widths. We see that the 
amplitude of the response typically increases at smaller 
scales indicating a change in the duration of the recom- 
bination epoch (i.e., the effective width of the visibility 
function). The Ci response also has an oscillatory com- 
ponent similar to a change of the position of the visibil- 
ity peak. These oscillations are most noticeable for the 
perturbations close to the visibility peak [z ~ 1100) and 
together with the overall slope of the Ci response demon- 
strate that the perturbation effectively changes the du- 
ration as well as the redshift of the recombination. 

2.3.2. Non-localized basis functions 

In this part we expand the perturbations in terms of 
two non-localized basis functions which, unlike the lo- 
calized case (with finite number of basis functions), do 
not suffer from non-orthogonality. If sufficiently many 
functions are taken into account, these basis functions, 
similar to the localized ones, can also be considered com- 
plete in practice. 

The most commonly used set of non-localized basis 
functions is the Fourier series: 



ipi{z) cx cos(«7r?;) 
^i{z) cx sin(i7r?/) 

lid 



i=0,l,2,. 
^= 1,2,... 



(8a) 
(8b) 



y 



Az/2 



where Az and Zmid are the width and central point in 
the redshift range of interest. Thus we have \y\ < 1 as is 
required for Fourier expansion. 

Alternatively, we can use Chebyshev polynomials of the 
first kind, T^, to form the basis. These modes are con- 
structed using the recursion formula: 

Ti+i{x) = 2xTi{x) - Ti^i{x) 

* http://camb.info/ 
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Fig. 2. — Perturbations in the Xc history, in the form of M4 splines (left) and the derivatives of the Cgs with respect to the amplitude 
of each perturbation (TT power spectrum in the center, EE on the right). 
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Fig. 3. — Similar to Fig. [2]but for perturbations in the form of Chebyshev polynomials. 



with initial conditions To{x) — 1 and Ti{x) - 
case, the variable x is replaced by y as given by Eq. (8c 
Chebyshev polynomials of the first kind are orthogonal 
with respect to the weight function w{x) = 1/ \/l — x"^ . 

The non-local basis functions are very different in na- 
ture from the localized functions discussed in § |2.3.1[ 
Therefore the response of the observables (here the C^'s) 
to the perturbation 5u{z) in the form of these functions 
is also expected to be rather different. Figure [3] shows 
the Ci responses when perturbing the ionizationhistory 
using Chebyshev polynomials with different frequencies. 
We see that perturbations with low frequencies, covering 
a large redshift range, lead to responses with much 
larger mean amplitudes when compared to the perturba- 
tions in the form of local bumps (Fig. [2]). However, as 
the frequency of the oscillations of the basis function in- 
creases, the response becomes weaker and its oscillations 
damp away. That is because neighbouring oscillations 
lead to similar responses in the C^'s with opposite signs 
and can partially cancel out each other. In other words, 
the CMB power spectra are less sensitive to the very high 
frequency perturbations in the ionization history. 

In principle, in the limit of large mode number, all 
bases work well. However, we have found that for the 
recombination history, although non-localized basis sets 
have their virtues, the z-localized bases are better, es- 
pecially if we are trying to describe narrow features in 



In this redshift. We will return to this point in § |3.1 



2.4. Constructing the eigenmodes 

We have so far introduced our choice of parametriza- 
tion for characterizing potential perturbations to the ion- 
ization fraction and illustrated how the perturbations in 
the form of these parameters affect the CMB power spec- 
tra when expanded in different bases. In principle all the 
qiS (in Eq. Q) are needed for a (nearly) complete re- 
construction of a general perturbation 5u{z). However, 
in practice data cannot constrain the perturbations in 
detail in many cases. As we saw in Fig. [3j very high 
frequency perturbations are expected to lead to much 
smaller signals. 

2.4.1. Eigenmodes with fixed standard parameters (XeMs) 

To avoid dealing with many correlated (and possibly 
weakly constrained) parameters (i.e., the (?i's), we con- 
struct a set of their linear combinations which are un- 
correlated with each other and only keep those combina- 
tions that are most constrained by data. This procedure 
provides a hierarchy of mode functions and their corre- 
sponding signals in the CMB temperature and polariza- 
tion power spectra. Exclusion of the weakly constrained 
eigenmodes does not affect the rest of the measurements 
since the eigenmodes describing the recombination per- 
turbations are by construction uncorrelated. 
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If N parameters are used to characterize the pertur- 
bations, i.e., 1 < i < N , the N uncorrelated parameters 
will be determined by the eigenmodes of the {N x N) 
Fisher information matrix: 



functions we started with: 



9^ Inpt 



where in the language of Bayesian analysis, pf = 
p(q|d, T) (with q — (gi, gjv)) describes the posterior 
probability of the parameters g^'s for the given data d in 
the theory space 7", i.e., an update from the prior proba- 
bility distribution of the parameters, pi — p(q|T) driven 
by the likelihood pf = £{q\d,T)pi/£ where C{q\d,T) = 
p((i|q, T) and the evidence £ = p{d\T) is a normalization 
factor. We include the T in the notation only if there is 
ambiguity in the theory space under consideration. 

Under the assumption of uniform priors for the qi's, 
the Fisher matrix reduces to: 



F- — — 



The derivatives are calculated at the fiducial values of 
the parameters, in this case gi = ... = (/at = 0. The 
ensemble average (..) is over realizations of the CMB 
sky and instrument noise. In the standard CMB anal- 
ysis with Gaussian signal and noise, we have C — 
exp(-Al'C"iA/2)/y27r|C|. Here A represents the tem- 
perature and polarization maps including both CMB sig- 
nal and instrumental noise and C — (AA^^) is the the- 
oretical pixel-pixel covariance matrix. With this likeli- 
hood function, the Fisher matrix simplifies to: 



F, 



2 



i9C 



In the limit of full sky observation, or in cut-sky cases 
where coupling between modes of different scales can be 
ignored, F can be written as: 



/sky 2^ — ^ — ?^ — ^ 



with 



£=2 



(9) 



e 

,22 



N 



where we have included CMB temperature T, i?-mode 
polarization and their cross correlation TE. Here nJ''^ 
stands for instrumental noise contribution to the power 
spectra and a is the width of the Gaussian beam. The 
effect of incomplete sky coverage has been naively taken 
into account by the /sky multiplier which reduces the 
effective number of observed modes . 

The Fisher matrix for N parameters, as any other 
N X N real symmetric matrix, has N independent eigen- 
vectors which can be chosen to be orthogonal to each 
other and normalized to one. So F can be decomposed 
as F = Sf S"^ where the columns of S are the eigenvectors 
of F with their corresponding (non- negative) eigenvalues 
on the diagonal of the real diagonal matrix f . The eigen- 
modes we are looking for can now be constructed us- 
ing these eigenvectors of the Fisher matrix and the basis 



N 



(10) 



If the (pi's happen to be orthonormal, then the eigen- 
modes Ek{z) will be: 



Ek{z)Ek'{z)w{z)dz = 6kk' 



(11) 



Here w{z) is the weight function with respect to which 
LPiS are orthonormal. Since Eq. (11) is not necessar- 



ily fulfilled, we enforce the Ek{zYs to be normalized to 
unity (as a matter of convenience) , which is equivalent to 
a renormalization of the eigenvectors of F. Although in 
general this could change the rank ordering of the modes, 
in our case a reordering was not required. Now, instead 
of the original fiS, the set of the eigenmodes can serve as 
basis functions for the expansion of perturbations (com- 
pare with Eq. (Ill 



5u{z) 



N 

E 

fc=i 



^lk Ek{z). 



(12) 



To be more specific, we will refer to these parameter 
eigenmodes which desc ribe perturbations to the his- 
tory by XeMs. In § |3.5| we will demonstrate two examples 
of perturbation reconstruction with diff erent numbers of 
eigenmodes taken into account (Fig. 10 ). We will see how 
well these eigenmodes serve as basis functions and also 
which features of the original perturbations are restored 
(or lost) if only a subset of the eigenmodes are used in 
the reconstruction process. 

The inverse square of the eigenvalues of the Fisher ma- 
trix can be used to forecast the error bars of the eigen- 
modes, i.e., fij = a^'^Sij, assuming the probability distri- 
bution of the parameters is multivariate Gaussian close 
to the maximum. For non-Gaussian likelihoods, the cr^'s 
give the lower bound for the errors. In the rest of this 
paper we use the term error for the Ci's, as the Gaus- 
sianity of the likelihood function close to its maximum is 
usually a good assumption. If the modes are sorted in de- 
scending order of eigenvalues, the first few (with smallest 
CTi's) will be most constrainable. Thus, the constrainable 
part of the perturbations to the ionization history can 
be described by the eigenmodes which have reasonably 
small uncertainties (i.e., high eigenvalues), while the rest 
is practically unconstrainable by the dataset under con- 
sideration. 

This procedure of using an orthogonal transformation 
to replace the parameters of the problem with a set of 
uncorrelated variables is called principal component anal- 
ysis or PCA for short. (The p arameter eigenmode s were 
used for CM B in Bond ( 1996 1 and subsequently in Bond 
et al. (19971 and many subsequent papers.) In general, 
the PCA needs to be applied to the whole history of 
ionization simultaneously (as well as the standard cos- 
mic parameters), since we do not know a priori how the 
ambiguity in one epoch affects the measurements of the 
perturbations in other epochs. However, if it turned out 
that a particular period of Xc history could be relatively 
well constrained, e.g. by other probes, one could leave 



8 



Far hang et al. 



that epoch out of the perturbations. Moreover, choos- 
ing a suitable parametrization, potentiahy changing over 
time to properly take into account the different physics 
at different epochs, is a necessary but not straightfor- 
ward task. In this work, the focus will be on the epoch 
of recombination since that is where the main CMB sig- 
nal is coming from. A more complete analysis for the 
whole ionization history or where different parts of it are 
considered simultaneously is for future work. 

2.4.2. Eigenmodes with varying standard parameters 

Above we assumed that the background cosmology was 
fixed and the parameter space only included perturba tion 
parameters, i.e., g^'s. However, as we will see in § |4.1[ 
the eigenmodes constructed this way do not necessarily 
stay uncorrelated with one another when standard cos- 
mological parameters are also being varied. Moreover, 
due to the correlation of the XeMs with one another as 
well as with the standard parameters, the constraints on 
the eigenmodes will possibly increase compared to the 
forecasted errors (reported in Tabled]). 

To avoid this issue, the eigenmodes of perturbations 
should be constructed in the presence of varying standard 
parameters. In this case the Fisher information matrix 
has the following general form 



where we have F 
andF, --/d'lnc 



F F 

^ ss ^ sp 



dsds 



(13) 



0^ ln£ 
dsdp 



■ pp — \ QpQp I , with s and p symbolically repre- 
senting the standard and perturbation parameters. To 
find the eigenmodes for perturbations after marginaliza- 
tion over the standard parameters, hereafter extended 
eigenmodes, or eXeMs, under the assumption of ellip- 
soidal parameter contours, we need to eigendecompose 
the perturbation block of inverse of the Fisher matrix, 
i.e., to find the eigenvectors of (F~^)pp. However, the in- 
version of the Fisher matrix is numerically problematic 
as it is ill-conditioned due the non-constrainable param- 
eters which exist in the parameter space under study. 

In the case of fixed standard parameters, the Fisher 
matrix reduces to the perturbation block only, and since 
the eigenvectors of a symmetric matrix and its inverse 
are the same (with inverse eigenvalues), there is no need 
to invert the Fisher matrix before its eigendecomposi- 
tion. Similarly to get the eXeM, i.e., the eigenvectors of 



(F ^^pp with F as in Eq. (13 1, we avoid the direct full 



inversion of F by noting that 

(F )pp — (Fpp — FpgFgj, Fsp) . (14) 

The eigendecomposition of Fpp — FpsFj/F^p then only 
requires the inversion of the well-behaved standard pa- 
rameter block. It is straightforward to directly check 
that (F~^)pp properly describes the marginal likelihood 
of the perturbation parameters: 

Here p and s are the arrays of the perturbation and stan- 
dard parameters. 



The eXeMs are then uncorrelated with each other (but 
not necessarily to the standard parameters) even in the 
presence of varying stand ard c osmic parameters. We will 
use this method later in § |3.8| to find the eXeMs describ- 
ing the perturbations to recombination history. Here- 
after, we use the term eigenmode with the very general 
meaning and reserve XeM only for the eigenmodes if they 
are constructed with fixed standard parameters. 

It is also possible to treat the standard parameters 
depending on the way they affect the power spectra. 
Among the standard parameters, yp has this unique 
property of influencing 's only through its impact on 
Xq. In other words, if we find the template in pa- 
rameter space corresponding to Fp, i.e., dXe/dYp, small 
changes in Yp can be mimicked by properly changing the 
amplitude of this template while Yp is left unchanged. 
Note that changes in other standard parameters either 
directly lead to changes in the CMB with no infiuence 
on Xc (such as rtg and ^s), or have both direct and 
indirect (i.e., through ATo) impacts on the CMB (such 
as r^b)- Therefore, when including a Yp-like parameter 
in the analysis, we can restrict our search for the per- 
turbation eigenmodes to the part of X^ space which is 
uncorrelated to the X^ template corresponding to this 
parameter as described above. In this way if perturba- 
tions to Ac are initially described by A'' parameters, we 
will in the end have the Yp-like parameter or its associ- 
ated template as one parameter and A^ — 1 eigenmodes 
which are uncorrelated to the Yp-like parameter. These 
A^ — 1 eigenmodes together with the Yp-like parameter 
fully parametrize the original A^-dimensional space of the 
perturbations. However, for the purpose of this work we 
did not further explore this possibility. 

3. PERTURBATION EIGENMODES FOR RECOMBINATION 

In this section we follow the procedure explained in 
§ |2.4| to find the eigenmodes for perturbations in the ion- 
ization fraction at high redshifts. For the most part we 
leave the standard parameters fixed, but return to the 
question of how these affect the mode functions later. 
We choose the redshift range of [200,3000] which cov- 
ers hydrogen and singly ionized helium recombination 
(z ^ 1100 and z ^ 1800 respectively) as well as part 
of the dark ages while being well above the standard 
reionization scenarios (z < 30). We assume the fiducial 
reco mbi nation history is given by the SRS, as explained 
unless otherwise stated. 



2.1 



In the following we compare the eigenmodes generated 
by using various bases and study some of the aspects as- 
sociated with them, such as their convergence and fidu- 
cial model dependence. Some consistency checks will also 
be presented. Special attention will be given to pertur- 
bation to helium recombination. We will also study how 
including the standard cosmic parameters in the analy- 
sis would change the eigenmodes of perturbations to the 
ionization history. At the end we propose an information- 
based criterion for cutting off the eigenmode hierarchy to 
be used in the data analysis. 

3.1. XeM construction using different bases 
In this section we take 6u{z) — (5 In Ac and try the 



five different sets of basis functions described in § 2.3 



Chebyshev polynomials and Fourier series as orthogona 
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Fig. 4. — The six most constrained XeMs for five different basis functions (with 160 parameters). The maximum and width of the 
Thomson visibiUty function have been marked in all figures. 



TABLE 1 

The forecasted standard deviations of the first six XeMs from the 
Fisher analysis for different observational cases. 



XeM 


1 


2 


3 


4 


5 


6 


CVL(£niax = 3500) 


0.003 


0.009 


0.013 


0.016 


0.022 


0.047 


CVL(£naax = 2000) 


0.011 


0.019 


0.024 


0.041 


0.094 


0.190 


CVL(^niax = 3500, T only) 


0.004 


0.021 


0.064 


0.103 


0.208 


0.275 


Planck-ACTPol{^max = 3500) 


0.015 


0.047 


0.068 


0.13 


0.22 


0.31 



non-local functions of redshift, and A/4 splines, triangu- 
lar and Gaussian bumps as localized functions. For the 
latter three, the width of the bumps is chosen to be inde- 
pendent of redshift. We choose o-j = Sz/2 for Gaussian 
and triangular bumps (Eqs. ([5]) and (|6|) and h = 1.5dz 
for M4 splines (Eq. Q). In all cases, dz = Az/{N + 1) 
is the spacing between the centers of adjacent bumps, 
where Az is the redshift range of interest and N is the 
number of basis functions used. For each set of basis 
functions we calculat e th e N x N Fisher information ma- 
trix as explained in § |2.4| where the N parameters are the 
amplitudes of the perturbations in the form of the basis 
functions, i.e., g^'s in Eq. (HI). The standard cosmic pa- 
rameters are fixed to their nducial values. For the data 
we simulate the T, E and TE spectra up to / = 3500 for 
a full-sky, cosmic variance-limited (hereafter CVL) CMB 
experiment, unless otherwise stated. We then construct 
the Fisher matrix (Eq. ([9|) and from it the N XeMs 
(Eq. ( 10 )). The first six XeMs for N — 160 are shown in 



Fig. |4] The first row in Table [ll shows the forecasted er- 
rors of these XeMs, obtained from the eigenvalues of the 
Fisher matrix. Note that including standard parameters 
in the analysis, e.g. MCMC simulations, can increase the 
error bars, as we will see later in § [4j 

We can see that the first six XeMs - which are the 
most constrained modes - all have the strongest varia- 
tions close to the maximum of the Thomson visibility 
function. The freeze-out tail is not perturbed signifi- 
cantly, and the oscillations around helium recombination 
(z ^ 1800) have much smaller amplitude than those at 
z ^ 1100. This is expected since the CMB anisotropies 
are most sensitive to perturbations during maximum visi- 
bility and features at low and high redshift are not weigh- 
ing as much in the CMB power spectra, once uncertain- 
ties close to z ~ 1100 are allowed. This in turn implies 
that only once the ionization history during hydrogen re- 
combination is known well can small modifications in the 
freeze-out tail or at z > 1400 be constrained. 
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We can observe another aspect of the XeMs. The larger 
the expected error bar the more small scale structure 
or high frequency oscillations the modes have and the 
further away from the visibility peak they probe. This is 
again understandable, since neighbouring ups and downs 
in the mode functions lead to partial cancelation of the 
effect on the Ci's. Once several oscillations are occurring 
close to z ~ 1100, signals produced farther away from 
maximum visibility can start competing with those from 
z ~ 1100, and hence become constrainable by the data. 

We also see that the first few XeMs are practically the 
same independent of the chosen expansion basis (Fig. [4]) , 
although individual perturbations in different bases lead 
to totally different C( responses (Figs. |2] and l3|. More- 
over, these modes are converged and do not change by 
including higher modes, as we will see in § |3.2[ However, 
as we go to modes with higher uncertainty (not shown 
here), the XeMs from different bases start to slightly dif- 
fer from each other. A larger number of basis functions 
are required to make these higher XeMs agree as well. 
Moreover, we found that the higher (poorly constrained) 
XeMs, in particular from the extended basis expansions 
such as Fourier, become dominated by numerical noise. 
The reason is that for weakly constrained modes where 
the higher frequencies start to play a more important 
role, the impact of adjacent ups and downs from the 
high frequency perturbations (e.g. sine functions) may 
not be well resolvable in the C^'s, resulting in their net 
effect being dominated by numerical noise. For the local- 
ized basis functions, as long as the individual bumps are 
numerically resolvable, we do not find this issue, because 
each perturbation has just one bump with no destructive 
neighbour. 

For more precise computations of the higher XeMs im- 
provements of the numerical treatment in Camb would 
become necessary. We tried several obvious modifica- 
tions, as well as different settings for the accuracy level, 
but were unable to stabilize the results for very high fre- 
quency modes. However, since in the analysis we are 
hardly using more than a few XeMs, for the purpose of 
this work this was sufficient. 

3.2. Convergence test 

In all cases considered above the number of modes was 
chosen to be iV = 160. We tested the convergence of the 
XeMs by trying N = 40, 80, 160 and 320 in different 
bases and found that by = 160 the first few modes 
are converged. Fig. [5] shows a representative example 
of this convergence with M4 splines as the basis. Thus 
the few constrainable XeMs are independent of the very 
high frequency (or highly localized) perturbations and 
including basis functions of higher order will not affect 
these eigenmodes. 

For the case of M4 spline functions the robustness of 
the results should also be checked against increasing the 
width of the kernel. By compari ng t he (first six) XeMs 



with h = l.bSz (as explained in § 3.1 1 to those with h = 
36z, we conclude that the modes have already converged 
for h = 1.5Sz and thus we will adopt h = 1.56 z in the 
rest of this paper. 

3.3. Fiducial model and dataset dependence 

It is important to note that the XeMs are by construc- 
tion fiducial model dependent. In principle, the observ- 



ables (such as C^'s) for different fiducial models respond 
differently to the same perturbations depending on the 
strength of the signals, at different redshifts, from the 
unperturbed fiducial model. 

As an example, in Fig. [6] we compare the eigenmodes 
for three fiducial histories. Two of the models have 
different CMB temperatures and in the third one lens- 
ing has not been included. In the first two, the different 
Tcmb's lead to different fiducial X^s. Here, the main dif- 
ference in the eigenmodes is their shift towards lower z's 
for the case with higher CMB temperature. This is con- 
sistent with the delayed recombination shown in Fig. [T] 
remembering that XeMs are primarily localized around 
the maximum of visibility where the C^'s are most sensi- 
tive to. For the latter case with no lensing, although X^ 
and the physics around recombination have not changed, 
there are still slight changes in some of the XeMs as seen 
in Fig. [6) 

We also checked the robustness of the eigenmodes 
against changes in the fiducial value of the parameters 
and the assumed reionization scenario. We tried a dif- 
ferent value for fib, as the parameter most strongly af- 
fecting the ionization fraction, \a away from its fiducial 
value. For the late reionization we tried an extended 
reionization scenario (i.e., X^ = 1 for z < 6; X^ = 0.5 
for 6 < z < 30 and X^ = Q elsewhere) radically different 
from our sharp fiducial reionization model (the default 
in Camb). For both of these tests the first six eigen- 
modes were found to be the same as our main eigenmodes 
(Fig 4| with tiny differences in the fifth and sixth modes 
for the latter case. 

This implies that, although the eigenmodes are fiducial 
model dependent, the constrainable ones are not practi- 
cally sensitive to changes in the fiducial model or its pa- 
rameters in the limits currently allowed by the data for 
the standard model of cosmology. That is because small 
changes in the fiducial parameters and the correspond- 
ing small changes in the ionization history only affect the 
XeMs at second order. Here by small we mean changes 
that lead to (smaller than or) the same order of magni- 
tude signal in the simulated data as the (few best) XeMs. 
The higher XeMs with larger uncertainties are more af- 
fected by the same changes in the fiducial parameters, 
as these changes are no longer considered small relative 
to these poorly constrained XeMs. This non-sensitivity 
of the best modes to the fiducial values of parameters 
does not contradict their significant correlation once the 
standard parameters are also allowed to vary, as we will 
see later in § |4.1[ 

We also studied the dependence of the XeMs on some 
properties of the simulated CMB datasets used for their 
construction, such as different Imax corresponding to the 
smallest scale information present in the data, and dif- 
ferent experimental noise levels. The results for a CVL 
experiment up to ^ = 2000 in temperature and polariza- 
tion and also a CVL experiment only sensitive to tem- 
perature (up to ^ = 3500) are shown in Fig. [t] As a more 
experimentally motivated case, we calculated the XeMs 
for simulated Planck-like datcPl (using 100, 143 and 217 
GHz channels, with effective galaxy-cut sky coverage of 
75%) and ACTPol-like data, including both wide and 

« http://www.rssd.osa.int/SA/PLANCK/docs/Bluebook-ESA- 
SCI282005291_V2.pdf 
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CMB temperature and in the other model lensing is not included. For the case of the two different CMB temperatures, the major difference 
is the shift in the eigenmodes associated to the shift in the fiducial and visibility functions (see Fig.Q. 

modes, determined from the eigenvalues of the Fisher 
matrix, become considerably larger. All these errors are 
calculated with the standard parameters fixed. However, 
the considered cases illustrate the general behaviour of 
the method. Taking into account the correlation between 
the perturbations and the standard cosmic parameters 
leads to relatively higher error bars, depending on the 
dataset used. We will return to this point later (table [s]). 



deep surveys ( Niemack et al.|[2010| ). As shown in Fig. [S] 
there is a tiny shift in the first mode relative to the mode 
for an ideal experiment and the changes grow as we pro- 
ceed to higher modes. 

More significant than the small changes in the XeMs 
constructed with different assumptions about data, are 
the forecasted error bars in different cases (see Table [1]). 
By removing the temperature at high £'s or the polariza- 
tion spectrum, the constraints on the amplitudes of the 
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3.4. Alternative parametrizations for perturbations 

We now turn to the parametrization defined by Eq. ([s]), 
Su{z) = S \n(X(. + a) with a > 0. As mentioned earher, 
this parametrization allows us to focus the perturbations 
on the freeze-out tail which are more difficult to recover 
with 6u{z) = Sln{Xc). This is because choosing cr > 
results in over-weighting the signal from the perturba- 
tions at low z (with low Xg) compared to a case with 
(7 = for the same value of Su{z). 

Fig. [9] shows the six most constrained XeMs with 
du{z) = dln{Xc + cr) for different values of cr, all 
constructed by using 160 Gaussian bumps as the 
basis functions and assuming CVL experiment up to 
^max = 3500. We see that the amplitude of the XeMs 
increases in the f reeze -out tail as a increases. It will 
be illustrated in § |3.5| that a relatively high value for cr 
is the preferred choice for studying perturbations that 
most significantly alter the freeze-out tail. 



3.5. Consistency check: eigenmodes as a complete basis 

Any function Xc{z) (in the redshift range under consid- 
eration) can be expanded in terms of these XeMs unless 
it has highly localized features compared to the highest 
frequency present in the basis functions or to the width of 
the bumps in the case of localized modes. That is because 
the XeMs are just linear combinations of the original ba- 



sis functions, and thus cannot have frequencies higher 
than the maximum frequency present in the basis func- 
tions. On the other hand, as is clear from Fig.[4j strongly 
localized features in possible perturbations to recombi- 
nation history are not constrained with CMB datasets. 
Therefore the lost features of an ionization model via 
expansion by these eigenmodes will not be measurable 
even if modes with higher frequencies are included in the 
analysis. In other words, the XeMs serve as a complete 
basis for the expansion of constrainable features in the 
possible perturbations in the recombination history. 

To demonstrate the reconstruction of perturbations us- 
ing the XeMs we choose two physically motivated ioniza- 
tion perturbations, one associated w ith physical correc- 
tions to the recombination process (Chluba & Thomas 
( |2011| , hereafter CT2010) and the other due to energy 
injection coming from a mod el of dark ma tter annihila- 
tion (using the description of ChlubapOlO ) . 

3.5.1. Standard recombination corrections 

The modification to Xc corresponding to CT2010 is 
shown in the top left panel of Fig. 10 (black solid line). 
This correction should be added to me Xc from the orig- 
inal version of Recfast (or the X^ from Recfast vl.4.2 
setting Heswitch = 0). At high redshifts one can see the 
effect of accelerated helium recombination caused by ab- 
sorption of photons in the Lyman continuum of hydro- 
gen. During hydrogen recombination the corrections are 
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Fig. 9. — The six most constrained XeMs with 5u{z) ■ 
bumps as the basis functions. 



: S ln{x + a) as the parameter, having three different values for a, and with Gaussian 



TABLE 2 

The projection of the modifications to recombination 

HISTORY ON THE FIRST SIX XeMS. 



XeM 


1 


2 


3 


4 


5 


6 


CT2010 
DM annihilation 


-0.32 
-0.31 


0.08 
-0.30 


0.16 
0.46 


0.02 
-0.14 


-0.09 
0.33 


0.25 
0.88 



caused by detailed radiative transfer effects as well as 
two-photon and Raman scattering events. The freeze- 
out tail is slightly higher than obtained with Recfast 
because of deviations from statistical equilibrium in the 
angular momentum sub-states. We note that with Rec- 
fast vl.5 a large part of all these corrections can be 
accounted for, however, these corrections are not explic- 
itly modelled using a physical description but have been 
fudged to reproduce the results obtained with detailed 
recombination codes. 

We project this Sln{Xc) on the 160 XeMs constructed 
from perturbations in the recombination history in the 
form of Gaussian bumps and with the perturbation pa- 
rameter being 6u{z) = S\n{X c) fo r the CVL case with 
^max = 3500, as described in § |3.1[ The figure compares 
the reconstructed perturbation for three cases with differ- 
ent number of XeMs included. First note that by includ- 
ing all 160 XeMs the original perturbation is practically 
fully recovered. If only the 15 most constrained modes 
are included, the hehum correction {z ^ 1800) and also 
hydrogen correction around z ~ 1100 are well restored 
while for lower z regions higher modes are required. The 



reconstruction by six XeMs, however, is most sensitive to 
variations around z ~ 1100 and cannot tell much about 
the helium correction. The projection coefficients for the 
first six XeMs are shown in Table [2] For this particular 
model of corrections to the perturbation scenario we see 
that the XeM 1, 3 and 6 are strongly dominant among 
the first six modes. 

The lower left panel in Fig. [lO] illustrates the relative 
difference between the temperature power spectrum for 
the reconstructed perturbations and the original full cor- 
rections. We see that with only six modes the error in the 
recovered C^'s is less than 0.1%. Remembering that the 
changes in the Ce's due to the full corrections are about 
a few percent, this shows that the main corrections to 
the CMB power spectra can be captured by just intro- 
ducing a small number of modes. The CMB data indeed 
are not sensitive to all the details in the freeze-out tail of 
recombination and during helium recombination, unless 
prior knowledge renders uncertainties at z ~ 1100 very 
small. As we will see below, part of the corrections from 
higher modes are compensated for by biasing the XeMs 
included in the analysis. 

3.5.2. Dark matter annihilation scenario 

As the second example we chose the perturbations aris- 
ing from a model of dark matter annihilation. It was 



computed using the description of Chluba (2010) with 
an annihilation efficiency /dm ~ 2 x 10~^^ eV/s. The 
difference with respect to Recfast is shown in the right 
panel of Fig. [T0| In contrast to the previous case, the 
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perturbations here are not concentrated around the max- 
imum of differential visibiUty but are most significant at 
lower redshifts. Therefore, for the decomposition of the 
DM p erturbations we choose Su{z) — 51n(Xe + 0.01) (see 
3.4 and Fig. |9]) to allow a better recovery of the rela- 



FlG. 10. — The reconstruction of two physically motivated Xc perturbation scenarios on (different number of) XeMs generated with 
Gaussian bumps (top) and the relative difference in the temperature power spectrum between the reconstructed perturbations and the full 
corrections (bottom). Right: The perturbations come from deviation from physical corrections to the recombination process (CT2010). 
Here the perturbation parameter is 5ln{Xc). Left: The perturbations are due to a model of dark matter annihilation. As the perturbation 
parameter we used <51n(Xc + 0.01) to better accommodate for the freeze-out perturbation. A case with 5u(z) = Sln{Xc), i.e., cr = 0, is 
shown for comparison. 

with Su{z) — 6 \n{Xc), i.e., with cr = 0, and with 15 XeMs 
taken into account is also shown. As expected, this recon- 
struction is much poorer compared to the previous case 
with Su{z) = S ln(Xo -1-0.01) due to its lack of coverage of 
corrections in the freeze-out tail. This demonstrates that 
when there is prior knowledge in favour of the freeze-out 
tail of recombination being affected, a parametrization 
with a > should be used in the analysis. However, 
it is still correct that the main signal is produced by the 
modifications close to z ~ 1100, even if the freeze-out tail 
apparently has the largest deviation from the SRS. This 
is why the first few mode functions for cr = do not have 
any strong low redshift tails. The eigenvectors naturally 
order the perturbations in the strength of the associated 
change in the CMB power spectra, as explained in § |3.1[ 
This point is visible from the lower right panel where the 
Ci difference is plotted for reconstructed perturbations 
with different number of modes included compared to the 



lively large perturbations in the freeze-out tail without 
the need to include too many modes. This procedure can 
be interpreted as placing a strong prior on (physically) 
expected changes in the freeze-out tail. 

The top right panel of Fig [TO] shows the reconstructed 
perturbation including three different number of XeMs. 
Here the recovered curve becomes very close to the orig- 
inal perturbation by including the first 15 XeMs, while 
six XeMs have a poor recovery of the low-z part. Note 
that the plots are illustrating 5\n{Xc{z)) although the 
XeMs and thus the decomposition of the perturbation 
are all performed with 6\n{Xc + 0.01). 

For comparison the reconstruction of the perturbation 
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Fig. 11. — The Svis = (vis — visflj) (left) and the relative changes in the TT and EE power spectra (middle and right) for the six most 
constrainable XeMs. 



full perturbations. Similar to the previous case, these 
differences are several times smaller than the changes in 
the Cf's caused by this model of DM annihilation, again 
meaning that these few modes can well capture the con- 
strainable features of the perturbations. 

Also if we look at the decomposition of the recombina- 
tion correction into the first six XeMs (see Tablej2]) we see 
that they all have comparable contributions. This seems 
reasonable if we remember that the mode functions, de- 
spite being weighted toward the low redshift part, still 
have a significant component at high redshift which need 
to be cancelled out to recover this pattern of perturba- 
tion with its low redshift modification. Therefore the 
neighbouring modes have the same order of magnitude 
amplitude to properly cancel out the high redshift per- 
turbations. This difference in the amplitude of the modes 
in principle allows us to distinguish this type of pertur- 
bation from those of CT2010. 

3.6. Impact of the eigenmodes on differential visibility 
and CMB power spectra 

It is worthwhile to see how the XeMs affect the visibil- 
ity function and the CMB power spectra. The left panel 
in Fig. 11 shows the change in the visibility function (nor- 
malized to the maximum of the fiducial visibility) for the 
first six XeMs. It is remarkable that relative changes in 
Xo of a few percent close to the maximum of visibil- 
ity, which lead to measurable effects in the CMB power 



spectra, only cause relative changes in the visibility of 
the order of 10"'*. This confirms the high sensitivity of 
the CiS to ti ny changes in the visibility. 

From Fig. 
mode, XeM 



11 we also see that the most constrained 
Thas an effect on the CMB power spectra 
consistent with changes the width of the visibility func- 
tion and a slight shift of its peak position. The second 
mode corresponds primarily to a shift in the visibility 
peak. The higher XeMs, with their several peaks and 
valleys, lead to less trivial changes in the width and po- 
sition of the visibility function. 

The relative changes in the CMB temperature and 
.E-mode polarization power spectra due to the first six 
XeMs are illustrated in the middle and right panels of 
Fig. [TT] To aid visual comparison, the amplitudes of the 
XcMsnave been chosen equal to their associated la's so 
that they would lead to comparable changes in the Cf's. 
It is interesting to note that for XeM 1, due to its nar- 
rower visibility width compared to the fiducial model, 
the high £ damping in the temperature and polariza- 
tion anisotropics is smaller. At the same time, because 
of fewer scattering opportunities for the photons, this 
mode leads to less polarization (negative SCf at low ^'s). 
The Ce signal associated with the second mode has the 
strongest oscillatory behaviour, consistent with a shift in 
the position of the visibility function. The mainly oscilla- 
tory impact of XeM 4 and 5 on the C^'s also suggests that 
an effective shift in the position of the visibility function 
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Fig. 12. — The six most constrained eXeMs. The solid blue lines correspond to modes constructed after marginalization over six standard 
parameters while for the dashed red curves Yp is marginalized over in addition. 

is present, while for mode XeM 3 and 6 the tail of the 
power spectrum is strongly damped (thus the negative 

SCj'^ at small scales), again corresponding to a change 
in the effective width of the visibility. 



3.7. Perturbations to helium recombination 

As mentioned early in this section, the redshift range 
chosen in our analysis of perturbations to ionization frac- 
tion includes the recombination of singly ionized helium. 
Some of the most constrained XeMs we found also extend 
up to z ^ 1600. These, therefore imply some impact from 
the helium recombination epoch on the XeMs. 

One way to confirm this statement is to limit the red- 
shift range of perturbations to mainly include singly ion- 
ized helium recombination, e.g., [1500,3000], while the 
total Xe (from both hydrogen and helium) is perturbed. 
We observe that the XeMs constructed this way have 
comparably large values at the lower redshift boundary 
{z = 1500) and would steeply go to zero if enforced by the 
imposed boundary conditions. This indicates that de- 
spite being restricted to the helium recombination epoch, 
the XeMs are still most sensitive to changes in the sig- 
nal from the hydrogen recombination and changes in X^. 
due to the helium recombination are hardly constrain- 
able, unless a properly chosen non-uniform prior on SX^. 
is imposed. 

As already emphasized, the parameters du{z) only 
characterize relative changes in X^ and the full descrip- 
tion of the ionization fraction depends also on the stan- 
dard cosmic parameters as well as the relevant the- 



oretical assumptions about the physics of recombina- 
tion. Among the standard parameters, Yp has a dis- 
tinct role in describing an aspect of the ionization frac- 
tion complementary to Su{z) by determining the maxi- 
mum total number of electrons available at each redshift: 

H 



.^c,max — -^cmax 



^^"max « (1 - Yp/2)Nk, where TVb is 



the baryon number density. Therefore, although Yp, in 
the first instance, requires to be marginalized over when 
constructing the eXeMs due to its intimate relation with 
the ionization fraction, it is also legitimate to treat re- 
combination perturbations and the maximum number of 
electr ons a vailable at each redshift separately. 

In § |4.1| we will use MCMC to measure constraints on 
Yp alongside the six standard parameters and the first 
few XeMs using (simulated) CMB data. Also in the next 
section, we explore how the eigenmodes change if they 
are marginalized over Yp. 

3.8. eXeM construction 

The XeMs we discussed so far were constructed with 
non-varying standard parameters and therefore can be 
considered as the limiting case of zero errors on the 
standard parameters. However, as mentioned earlier in 
§ |2.4.2[ the eigenmodes become correlated when they 
are simultaneously being measured with the standard 
parameters, due to their degeneracy with the standard 
parameters. The strength of the impact of these corre- 
lations on the XeM estimation depends on the (prior) 
constraints on the standard parameters. It is therefore 
worthwhile to see how the modes and their rank ordering 
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Fig. 13. — The first three eXeMs with only Yp being marginaUzed over (dashed blue curves). For comparison, the first three XeMs are 
also plotted (solid black curves). 



TABLE 3 

The forecasted standard deviations of the first six eXeMs from the 
Fisher analysis constructed by marginalization over different number 
of standard cosmic parameters and for different observational cases. 

In all oases, £max = 3500. 



cXcM 


1 


2 


3 


4 


5 


6 


CVL, marg: six std 


0.011 


0.012 


0.029 


0.052 


0.059 


0.064 


CVL, marg: six std + Yp 


0.011 


0.027 


0.029 


0.052 


0.059 


0.071 


Planck-ACTPol, marg: six std 


0.058 


0.074 


0.189 


0.308 


0.439 


0.532 


CVL, marg: ris. As 


0.009 


0.011 


0.016 


0.018 


0.033 


0.059 



change if the standard parameters are allowed to vary as 
well. 



Figure 12 illustrates the first six eXeMs constructed af- 
ter marginalization over the main six and seven (includ- 
i ng Yp) standard parameters, constructed as described in 
§ |2.4.2| The eXeMs have stronger high redshift features 
compared to the XeMs. This implies that the degeneracy 
between the standard parameters and some features in 
the perturbations of the ionization fraction has pushed 
back some patterns of high significance to lower levels 
opening up the room for some high redshift or higher 
frequency patterns which only had the chance to show 
up at lower significance XeMs. 

The modes in the shown two cases (i.e., marginal- 
ized over six and seven standard parameters) differ only 
slightly. That is because Yp is rather weakly constrained 
using CMB data alone and in the presence of other stan- 
dard parameters its role in shaping the eigenmodes is 
only secondary. If, on the other hand, we hold the six 
standard parameters fixed and only let Yp vary, the eigen- 
modes will be more significantly affected (see Fig 13) 



The reason is that Yp is comparable in significance to 
small changes in the ionization fraction and marginaliz- 
ing over it, without the dominance of the standard pa- 
rameters, would lead to marked changes in the eigen- 
modes. The forecasted errors on the first six eXeMs with 
and without Yp included is compared in Table [3] We see 
that the errors are mostly the same, again implying the 
subdominant role of Yp. In terms of the errors the most 
affected modes are eXeM 2 and 6. 

3.8.1. Projecting eXeMs onto XeMs 



It is instructive to see how the eXeMs can be con- 
structed from the XeMs. Table H] shows the coefficients 
of the projection of the first six eXeMs (Fig. 12 1 on the 
best eight XeMs (Fig.|4]). Note that the most constrained 
eXeMs have their strongest projections along these first 
few XeMs and the contribution from all other modes is 
at most about a percent for these first six eXeMs. This 
means that allowing the standard parameters to vary 
mixes and rearranges the first few XeMs with negligible 
leakage from higher neglected XeMs. 

More specifically, the first eXeM is very similar to XeM 
2, with some contribution from XeM 3 and XeM 4. The 
second eXeM is mainly a combination of the first, third 
and fourth XeMs, and eXeM 3 and eXeM 4 have com- 
parable contributions from the first six and eight XeMs, 
respectively. The eXeM 5 and 6 are dominated by higher 
modes (XeM 6-8). As a result the space covered by the 
first few XeMs is not practically hugely different from the 
space covered by the first few eXeMs if proper number 
of modes (here eight XeMs for six eXeMs) are included 
in the analysis. This implies that using similar number 
of XeMs and eXeMs in an analysis of possible recombi- 
nation perturbations should give similar results for the 
reconstructed modification in the ionization history, at 
least for the CVL case where relatively large number of 
eigenmodes are included. The main advantage of the eX- 
eMs is that one can obtain more realistic estimates for the 
error bars directly after the construction of these modes. 
Also the eigenmode measurement is robust against in- 
cluding new eigenmodes in the analysis if the eXeMs are 
used as the modes. However, it also turns out that the 
eXeMs perform better than the XeMs in the simulated 
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TABLE 4 

The coefficients of projection of the six most constrained eXeMs on the 

FIRST EIGHT XeMS. 





XeM 1 


XeM 2 


XeM 3 


XeM 4 


XeM 5 


XeM 6 


XeM 7 


XeM 8 


eXeM 1 


-0.00 


-0.90 


0.21 


-0.29 


-0.12 


0.02 


0.02 


0.01 


eXeM 2 


-0.76 


-0.05 


-0.48 


-0.35 


0.16 


-0.02 


-0.02 


-0.01 


eXeM 3 


-0.31 


0.34 


0.42 


-0.38 


-0.54 


-0.41 


0.02 


-0.11 


eXeM 4 


-0.35 


0.13 


-0.48 


-0.31 


-0.50 


0.27 


-0.28 


0.37 


eXeM 5 


0.24 


-0.01 


-0.15 


-0.20 


0.06 


-0.48 


0.63 


0.49 


eXeM 6 


-0.19 


-0.22 


-0.07 


0.41 


0.00 


-0.41 


-0.55 


0.49 
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Fig. 14. — The best three eXeMs with the two inflationary parameters Ag and rig being marginalized over, compared to the XeMs, where 
all standard cosmic parameters are held fixed in the construction process of perturbation eigenmodes. 

analysis of Planck data, where only 1-3 modes seem to and supernova surveys (e.g. LSSTp^ Pan-STARR^^ 
be constrainable. BigBOSS, WFIRST{^ EUCLII^, the corresponding 

eigenmodes would be constructea after marginalization 
only over these inflationary parameters. 

Fig [T4| compares the first three XeMs with eigenmodes 
marginalized over and rig- The forecasted errors on 
these eigenmodes (from the Fisher analysis) are pre- 
sented in the last row of Table |3] Not surprisingly, these 
modes have smaller errors compared to the eXeMs which 
have been made after marginalization over six standard 
parameters, and have larger errors compared to XeMs 
(with no standard parameter varying) . These modes and 
similar ones after marginalization over different sets of 
standard parameters smoothly bridge the gap between 
the XeMs and the eXeMs. Depending on the datasets 
available at the time of real data analysis, the proper 
eigenmodes marginalized over the appropriate standard 
parameters must be constructed. With the current (and 
very near future) surveys, the most realistic choice are 
the eXeMs, constructed according to the experiment un- 
der consideration, which should be quite similar to the 
Planck-ACTPol-like case studied here. 



3.9. From eXeMs to XeMs 

The two sets of ionization perturbation eigenmodes in- 
troduced and constructed so far, i.e., the XeMs and the 
eXeMs, allow us to best describe and measure the uncer- 
tainties in the ionization fraction in the two extreme ends 
of our knowledge of the standard cosmic parameters. The 
eXeMs present a case where the tightest constraints on 
the standard (six) parameters are from CMB data alone. 
Therefore a simultaneous measurement of the standard 
parameters and the uncertainties in the ionization frac- 
tion, using the CMB dataset at hand, is required. The 
construction of the XeMs, on the other hand, assumes 
the standard cosmic parameters are measured with high 
accuracy from other cosmological probes and the CMB 
data are only used for the direct measurement of the 
ionization history. In other words, the XeMs, by ignor- 
ing the uncertainties in the standard parameters, extract 
the maximal amount of information that the CMB data 
would ever have to offer about the ionization fraction. 

Between these two limiting cases, there is a gray re- 
gion where, depending on the dataset at hand, tight pri- 
ors from non-CMB surveys can be imposed on some of 
the standard parameters while the rest are marginalized 
over when constructing the eigenmodes. For example, if 
all standard parameters, but the inflationary ones and 
rig, are measured to very high precision by other probes, 
such as large scale structure, baryonic oscillation, lensing 



3.10. Criteria for truncating the eigenmode hierarchy 

For the full reconstruction of perturbations to the ion- 
ization fraction, all eigenmodes are needed in principle 

http://www.lsst.org/lsst/ 

http:/ /pan-starrs. ifa.hawaii.edu/public/home. html 
http://wfirst.gsfc.nasa.gov/ 

http:/ /sci.esa.int/science-e/www/object/index.cfm? 
fobjectid=42266 
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since they form a complete basis set. In practice, se- 
quentially adding modes rank ordered in the (possibly 
renormalized) eigenvalues fk — cr^'^ of F from high to 
low gives a rapidly diminishing return once one goes be- 
yond a dozen or so. And often we can learn much from 
using just the first few. As more modes are added, the 
width covered by th e all owed trajectories increases, 
as Figs. [22] to [26] in § |4.3| show. The errors in those stan- 
dard cosmic parameters which are correlated with the 
X(. eigen-parameters also increase. On both counts, it 
behooves us to develop criteria for selecting which modes 
to keep, bearing in mind Occam's Razor for minimizing 
the number of new parameters to be added. Thus we 
show in § H] what happens when one mode, a few modes 
and a handful of modes are added. To be more quanti- 
tative, we explore a criterion based on not allowing the 
Shannon entropy to increase too much as the next eigen- 
modes in the hierarchy are added. 

The information action is defined in terms of the a 
posteriori probability of the variables pt and th e evid ence 
£ as 5i,f(q) = Inp^^ — ln£. Recalling from 



2.4.1 



the 

a posteriori probability pf = p{q\d,T) of variables q = 
((71, gAf) given the theory space T and the datasets d 
is related to the a priori probability Pi = p(q|T), the 
likelihood C{q\d,T) = p{d\q,T) and the evidence £ = 
p{d\T) through Bayes theorem: pi — C{q\d,T)pi/£ ■ The 
information action can then be written in terms of pi and 
C: 



5Lf(q) =lnpri +ln£- 



(15) 



For basic information th eoretic and Baye sian notions and 
notations see, e.g., the MacKay (20031 textbook. The 
framework given he re was used in a CMB context by 



Farhang et al. (2011 ).) For us the g^'s are the amplitudes 
ot the ordered eigenmodes for XeMs or eXeMs. Gener- 
ally the fluctuations in the standard cosmic parameters 
from their maximum likelihood values are included along 
with these eigen-parameters. We shall assume the prior 
distribution of the parameters to be uniform in the qk- 
The expansion of 6*1. f to quadratic order is the basic per- 
turbative approach used throughout this paper, leading 
to a Gaussian pf. 5i^f (q) ss iSi^m + q^rq/2 in terms of 
the Fisher matrix and the information action minimum 

Sl,m = -ln(pi>Cmax)- 

The posterior Shannon entropy is related to the final- 
state ensemble-average of the information action and the 
evidence: 

5f=(lnpf-i)f = (5Lf(q))f + lnf . 
For the quadratic order expansion it is 

Sf « ^TrlnF-i + ^ ln(27r) + ^Tr ((qq^)F) 



lTrlnF-i + ^(ln(27r) 



1) 



The second line follows from the first since the correlation 
matrix of the q is (qq^) — F~^. The associated evidence 
involves the information action minimum, ln£ « S{ — 

C N 
>->I,ra 2"" 

The entropy associated with mode n is 
^„ = In /„ -f (1 -f In 27r)/2 = S{< n) - 5(< n - 1) . 



It is a finite difference of the total entropy of the first n 
modes in the eigen-hierarchy, 

5(<n) = ^(l + ln2^)--^ln/fc 



fc=i 



and 



(s)„ ^ S{< n)/n 



is the associated mean entropy-per-mode. Fig. 15 shows 
how the relative entropy Sn — Si and the mean entropy 
(s)„ — Si = S{< n)/n — Si grow with n for the modes 
derived from the localized Gaussian expansion. We also 
plot two versions of " white-noise" entropy. 

5(wn, < n){a^) = n(lncr + (1 + ln27r)/2). 



variance cr 



k<r 



mean — weight — /fe /'^l ^ ■ 



k<.n 

These are entropies maximized subject to the constraint 
that we only have knowledge of the integrated cr^, 
whereas S{< n) is the maximized entropy given knowl- 
edge of the full spectrum {ff^^}- The mean- variance 
white noise lies above S{< n) and the mean-weight white 
noise lies below. The mean-weight behaviour is domi- 
nated by a ln(n) rise, since the total weight of modes 
below n, ^T^J2k<n fk' quickly approaches a constant, re- 
flecting the dominance of the high-weight eigenmodes in 
the sum. 

We first discuss why we do not use the traditional ev- 
idence ratio often used in Bayesian theory to decide if 
a new parameter qn should be added. The log-evidence 
difference for the addition of q„ is 

Aln£„ EE ln£(< n) - \n£{< n - 1) 

= 5„ - 1/2 - A5i,„,. 

This requires evaluation of the change in the infor- 
mation minimum. It also has the usual disadvan- 
tage of depending upon the fk measure. Although 
using eigen-parameters ensures the same dimensional- 
ity for the different /fc, it does not fully remove this 
re-parameterization ambiguity since there can be a k- 
dependent scaling. (In fact, we have usually renormal- 
ized our fk so that the associated eigenmodes Ek{z) have 
unit norm upon z-integration.) 

Our preferred approach for hierarchy truncation is to 
use suitably-defined entropy differences. In particular we 
wish to set a threshold control on the injection entropy, 

'5*S'inj,n = Sn , 

the entropy from adding mode n relative to the mean- 
entropy from all < n modes. It is related to the rela- 
tive increase in phase space volume V{< n) = exp(S'(< 
n) — n/2) = exp(n((s)„ — 1/2)) associated with mode 
additions: 



In 



F(<n + l)/F(<n)("+i)/" 



Sn 



We chose S'inj.n instead because it is zeroed out for mode 
one, but S'inj.n quickly approaches S'„+i — (s) . For exam- 



ple, if we impose a A^t ~ 1/2 threshold in Fig. 15 on the 
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CVL XeM case, we would use only one mode, whereas 
ASt ~ 1 picks up about 5, ~ 3/2 harvests about 10, and 
2 gives about a dozen. Similar tales can be told for the 
eXeM CVL case and for both Planck+ACTpol forecasts. 

Another more erratic measure is relative injection 
jumps, which is nearly Sn+i — Sn. In Fig. 15 the nega- 



10 



tive of this is plotted for clarity of presentation. Either 
reading off from the figure, or the using the lists of errors 
in Tables [l] and [sl the sample threshold A5t = 1/2 again 
yields only a mode or two. 

The fluctuating nature of Sn+i — Sn implies we can 
use it to split the modes into groups of similar informa- 
tion content which arise by thresholding it. If a mode 
is selected to be included in the analysis by, say, sharp- 
thresholding the injection entropy, it is logical that all 
of its co-modes be included, which is akin to softening 
the threshold. The groupings found with ASt = 1/2 im- 
posed upon Sn+i — Sn creates boundaries at one mode, 
five modes, and so on. These are, not surprisingly, simi- 
lar to mode numbers obtained as we move the threshold 
on injection up, hence that criterion can be used instead 
to define mode groups. 

Although these entropy difference criteria imply that 
relatively little additional information is gained by in- 
cluding more than a handful of higher modes, in real 
data analysis the situation is subtler, with other crite- 
ria important to consider. For example, depending on 
how close the assumed model is to the true underlying 
history, our measurements of standard cosmic parame- 
ters might be biased. In that case one would like to add 
enough modes to remove the bias, sequentially checking 
if the recovered values of the standard parameters are 
robust against the introduction of the next eigenmode. 
A reasonable strategy is to add one mode-group at a time 
to the analysis until the biases are removed. In the next 
section, we show how varying the mode number cutoff 
affects our results, roughly following this grouping pro- 
cedure. 

4. MEASURING THE AMPLITUDES OF PERTURBATION 
EIGENMODES FOR SIMULATED DATA 

Having constructed the eigenmodes, their amplitudes 
can now be considered as additional parameters to be 
plugged into CoSMOMCp] In this section we investigate 
how well the most constrained XeMs and eXeMs can be 
measured by simulated data. To study the impact of 
these new variables on the standard parameter estima- 
tion, we first consider the case in which th e dat a are both 
simulated and analyzed using the SRS (§ 2.1). We then 
study the case for which the effects of physical correc- 
tions to the recombination history (CT2010, see § 3.5 1 
are included in the constructions of the mock data, but 
are neglected in the fiducial recombination model used 
in the analysis. Here the question is how well the eigen- 
modes compensate for the deviations from the fiducial 
model and how much the data are telling us about the 
amplitudes of the modes. We also briefly discuss how the 
eigenmodes should be used in a more general case where 
little prior knowledge about the recombination perturba- 
tions is available. 

4.1. Case 1: The standard recombination scenario 
■'^ http:/ /cosmologist.info/cosmomc/ 
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Fig. 15. — Various measures of entropy differences defined in the 
text are plotted against increasing eigenmode number, for (a) the 
XeM case and (b) the eXeM case. The Cosmic Variance Limited 
mode results have heavy lines or points, the Planck+ACTpol fore- 
cast has lighter lines and points, as indicated. They look quite 
similar. For this figure, the modes are determined by the densely- 
packed Gaussian bump expansion, but the triangular and spline ex- 
pansions look similar, differences becoming notable only at higher 
n. The basic information on growing entropy is given by Sn — Si, 
and the mean difference S'<„/n — S\, with the latter bracketed by 
the two white-noise entropies. The criteria for threshold selection 
discussed in the text involve the injection entropy, Sn — S<n/n and 
Sn+i — Sn (plotted with opposite sign for clarity). Thresholding at 
ASt ~ 1/2 selects the first mode or two, but mode-groups are also 
evident, suggesting modes should be added in blocks rather than 
singly as we eliminate bias, check convergence and demonstrate 
robustness. 

As the first example we choose the fiducial recombina- 
tion model (here the SRS) to be identical to the ioniza- 
tion history used in the simulation of the data. We ran 
COSMOMC to estimate the best-fit values and errors of 
the six standard parameters together with those for the 
perturbation eigenmodes. We tried the two sets of eigen- 
modes described before: the first five XeMs and the first 
six eXcMs. The number of modes in each case was cho- 
sen in rou gh agreement with the mode cutoffs described 
The simulations were carried out for a CVL 



3.10 



experiment. 
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Fig. 16. — 2D contours for the amplitudes of some of the best five 
XeMs as measured by a CVL experiment and with the standard 
recombination history. Here the six standard cosmic parameters 
were also allowed to vary. That is why the measured eigenmodes 
are correlated. 




-0.040.02 0.02 





-0.0+0.02 0.02 



-0.1 0.1 




-0.0+0.02 0.02 



-0.1 0.1 




71 
70.8 
70.6 
70.4 
70.2 

70 




■6std 

■ 6std+ 5XeM 

■ 6std+ 6eXeM 



0.11 



Q._, h 

dm 



(^111 



0.112 



0.97 



0.965 



0.96 



0.955 




— 6std 

■■6std+5XeM , 

- 6std+6eXeM) 



0.0225 



0.0226 

b 



0.0227 



3.21 



3.2 



< 3.19 



o 



3.18 



3.17 




■6std 

■ 6std+ 5XeM 

■ 6std+ 6eXeM 



0.08 



0.085 0.09 



Fig. 17. — 2D contours for the amplitudes of some of the best six 
eXeMs as measured by a CVL experiment and with the standard 
recombination history. The six standard cosmic parameters were 
also allowed to vary in the simulations. 



Fig. 18. — The contours for the standard parameters as measured 
by an ideal experiment in the presence of five (six) XeMs (eXeMs) 
compared to the case with no eigenmodes included. The input 
value of the parameters is shown by the black diamond. 



TABLE 5 

The standard deviations of the first five 
XeMs from chains produced by CosmoMC, 
marginalized over the main six standard 
parameters. 



XeM 


1 


2 


3 


4 


5 


la 


0.046 


0.030 


0.057 


0.088 


0.086 



One expects no detection of eigenmodes since the fidu- 
cial model for Xa and the underlying model used to sim- 
ulate data are the same, as verified by Figs. [16] and [TT} 
Also, by construction there is almost no visible correla- 
tion betvifcen the measured parameters for the eXeMs, 
at least sufficiently close to the best-fit model where the 
assumptions of the Gaussianity for the likelihood surface 
approximately holds. However, Fig.[T6]indicates that the 
XeMs become partially correlated with each other, al- 
though by construction these were initially uncorrelated 
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Fig. 19. — The derivatives of the C„ ' 's with respect to some of the standard parameters. 



BThe reason is that the standard parameters were 
Id fixed during the process of XeM construction, but 
now that they are allowed to vary, their degeneracy with 
the XeMs induces correlations. These new correlations 
lead to larger errors than those deduced from the simple 
Fisher analysis ( Table [T]cf. Table [S]) and can also change 
the rank ordering of the modes, e.g., the error on XeM 2 
is smaller than XeM 1 (Table [5]). 

The standard parameters remain unbiased, as the 
model used for simulating data and the theoretical model 
used in the analysis were the same. This is no longer true 
once rec ombination corrections to the SRS are added (see 
Fig. 20). However, the correlations of the eigenmodes 



with some of the standard parameters incre ase the er 



we see 
ris and 



rors of the standard parametersj^ From Fig. [18 
that among the standard parameters, flhh^ , 
are the ones most affected by the introduction of the 
eigenmodes into the analysis. This can be understood 
by noting the relatively high degeneracy between these 
parameters and some of the eigenmodes. The most evi- 
dent one is the correlation of Ug with the first XeM which 
by changing the width of the visibility func tion lea ds t o 
a tilt in the power spectra (compare Figs 19 and 11). 
For the case of ^lt,h^ and As it is harder to give a visual 
interpretation, flhh^, leading to both tilt changes and os- 
cillations in the Ce's, correlates with most of the first five 
XeMs (the highest being with XeM 1), while As, being 
an amplitude multiplier, mainly correlates with XeM 1. 
These correlations between the standard parameters and 
the eigenmodes emphasize the fact that uncertainties in 
the recombination scenario in particular undermine our 
abili ty to measure the prec ise values of Ug and il^h^ (see 



e.g., |Shaw Ez Chluba 2011). Also note that the changes 
in the error bars ot the standard parameters are actually 
practically ind ependent of which set of eigenmodes are 
used (Fig. 18). This suggests that in terms of the stan- 



dard parameter estimation, the use of XeMs or eXeMs 

We confirmed this statement by running MCMC with non- 
varying standard parameters. 

^® It should also be noted that the correlations between the 
standard parameters themselves also change when the eigenmodes 
are introduced. 



should not lead to vastly different results in the parame- 
ter estimation. However, the perturbations are measured 
to higher accuracy with the eXeMs (Table [3| than XeMs 
(Table [5| especially if only a few modes are included 
in the analysis. Therefore, as long as only CMB data 
are used, the eXeMs are the more appropriate choice of 
eigenmodes. 

Finally, we studied how much the presence of pertur- 
bations to recombination could affect our ability to de- 
termine the precise value of Yp. The abundance of he- 
lium affects the CMB anisotropics mainly because more 
helium implies fewer free electrons during hydrogen re- 
combination. Consequently, Yp should also couple sig- 
nificantly to the perturbation eigenmodes. We therefore 
performed simulations in which Yp was also allowed to 
vary. The analysis was performed with three and five 
XeMs in the Planck-ACTPol-like and CVL case, and 
with six eXeMs for the simulated CVL data. Table |6] 
compares the la error bars on Yp in these cases. We see 
that for the CVL case similar number of XeMs and eX- 
eMs used as the eigenmodes lead to similar constraints 
on helium abundance. However, a Planck-ACTPol-like 
observation gives a few times larger error due to lack of 
very high sensitivity to very small scales, although fewer 
XeMs compared to the CVL case have been used. 

4.2. Case 2: A perturbed recombination scenario 

As the second example of parameter estimation and 
perturbation reconstruction, we simulate dat a ass uming 
the recombination calculation of CT2010 (Fig. 10), while 
we take the fiducial model to be as of Recfast vl.4.1 or 
older (equivalent to SRS with Hcgwitch — 0, to remove the 
helium correction which has been assumed as part of the 
perturbations in the data). The purpose here is to find 
out how well the biases in the standard parameters due to 
this lack of knowledge about the physical corrections can 
be removed by including the perturbation eigenmodes, 
and whether or not, data can reconstruct part of the 
true recombination history. 

Constructed from COSMOMC chains for a CVL exper- 
iment. Fig. 20 illustrates the 2D-contours of some of the 
standard parameters. The large biases in the estimated 
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TABLE 6 

Yp AND ITS MEASURED ERROR FROM SIMULATIONS FOR A CVL AND A PlANCK-ACTPOL-LIKE EXPERIMENT, WITH XeMs AND 
eXeMS TAKEN INTO ACCOUNT AS PERTURBATION EIGENMODES. COMPARED TO THE CASE WITH NO EIGENMODES. 



CVL(std) 


Planck-ACTPol(std) 


CVL(std + 5 XeMs) 


CVL(std + 6 eXeMs) 


Planck-ACTPol(std + 3 XeMs) 


Yp 0.240 ±0.0016 


0.240 ± 0.006 


0.239 ±0.005 


0.240 ± 0.004 


0.238 ±0.017 



0.97 
0.965 

0.96 
0.955 

0.95 
0.945 

0.94 





6stcl(CVL) 

- - - 6 std+ 8XeM (CVL) 
6 £tcl+ 20eXeM (CVL), 



0.0223 0.0224 0.0225 0.0226 0.0227 0.0228 



3.3 



3.25 



o 
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Fig. 20. — Contours of some of standard parameters for CT2010 case, with eight XeMs in one case and 20 eXeMs in the other case 
included in the analysis, compared to a case where no perturbation eigenmodes (of any kind) has been included (the solid red curves). The 
simulations are performed for a CVL experiment. The input value of the parameters is shown by the black diamond. 




Fig. 21. — Similar to Fig. |20| but for a Planck- ACTPol-like experiment. Here three eigenmodes were added for both the XeM and eXeM 
case. 



values of the parameters when only the six standard pa- 
rameters are measured is due to the mismatch between 
the ionization histiry in the theoretical model and the 
data. Here only contours for par ameters with the la rgest 
biases are shown. See also Shaw & Chluba (j2011|). To 
compensate for this mismatch we separately add to the 
parameters the two different sets of the eigenmodes, the 
XeMs and eXeMs, as the new parameters. 

As Fig. [20] demonstrates, this eliminates the biases 
in the stanHard parameters, however, at the cost of in- 
creased error bars. In particular for Vl\yh? the difference 
is large. The error (corresponding to ^ 5(t without eigen- 
modes) increases by a factor of ^ 5 when the (e)XeMs 



are included. For rig, similarly, the error (of ^ 7a) de- 
creases by a factor of 2. Also for we find a similar 
degradation, while for r the difference is rather small. 

Our computations indicate that with the XeMs as the 
eigenmodes and for a CVL observation, the minimum 
number of modes required to remove the bias from these 
standard parameters is six. However, we included eight 
modes in the analysis to ta ke into account the mode- 
selection criterion of § 3.10 (determined by the relative 



injection jumps). 

We also observe that the recovered values for the am- 
plitude of the XeMs are biased (compared to the theo- 
retically expected values from direct projection on the 
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XeMs, Table [2 1, and change by varying the number of 
modes included in the analysis. That is due to the cor- 
relation of the XeMs in the presence of the standard pa- 
rameters, and the fact that not all XeMs are included 
into the parameter estimation. As a results, parts of the 
perturbation that project on the neglected higher XeMs 
leak into the lower XeMs. The bias in the measured XeM 
amplitudes is similar to the bias in the standard param- 
eters when there are no eigenmodes in the analysis, but 
with a much lower significance. 

For the same reason, the errors on the XeMs also 
change when the number of modes included in the anal- 
ysis changes. However, as mentioned before, due to the 
low significance of the perturbation detection for most 
of the XeMs this is not as important as for the main 
cosmic parameters. For the CVL simulations with six 
and eight XeMs included, we see that the most signifi- 
cant contribution comes from the first mode (respectively 
fii ^ -0.23 ± 0.05, /ii = -0.18 ± 0.04) while the other 
modes are consistent with zero. This is also true for a 
Planck-ACTPol-like case, which we will come to shortly, 
where the first XeM is measured to be /ii = —0.22 ±0.06 
Hi = — 0.24±0.12 for one and three XeM measurements. 

If instead eXeMs are used as the perturbation eigen- 
modes, our computations show that at least ten modes 
should be added to get rid of the bias for a CVL case. 
However, as a test case, we tried including the best 20 
eXeMs (see Fig. [201). 

We also found that although the errors on the stan- 
dard parameters keep increasing by adding more eXeMs 
to the analysis up to around the tenth mode (which is 
required to remove the bias) it stays more or less the 
same afterwards. This suggests that in terms of the con- 
straints on the standard parameters, we do not lose much 
by increasing the number of eXcMs. Besides, including 
more eXeMs does not affect the measurement of the pre- 
viously included eigenmodes, as they are by construction 
uncorrelated (in the presence of standard parameters). 
Including more eXeMs, on the other hand, makes the re- 
constructed perturbation s cl oser to the input model of 
perturbations (as in Fig. 10 1. However, as the errors of 
modes increase by going to higher orders, the error on 
the reconstructed curve increases. We will address this 
point in the next section. 

Among the first 20 eXeMs for a CVL experiment, 
the modes with the most significant contributions (i.e., 
with at least Icr detection) are /i2 — 0.11 ± 0.02, fi^ — 
0.10 ± 0.03, HQ = -0.31 ± 0.16 and /in = -0.36 ± 0.24 
(compare to their theoretical prediction from direct pro- 
jection of the perturbations on the eXeMs: fi2 = 0.14, 
yU3 = 0.10, fi9 = —0.33 and fin = —0.39). The reason 
that the recovered value, though close, is not exactly the 
same as the forecast is that the assumption of the Gaus- 
sianity of the distributions of the eXeMs and the stan- 
dard parameters is only approximate. Also the eigen- 
modes have been slightly smoothed in the construction 
process which may cause numerical inaccuracy and in- 
duce slight correlation between the smoothed modes. By 
comparing the theoretical values of projection of the per- 
turbation on the eXeMs and their forecast errors (from 
Fisher analysis) we do not expect any perturbation de- 
tection after eXeM 11. 

Fig. [21] shows similar contours but for a simulated 
Planck^TSCTFol-like observation. For the analysis we 



used the eigenmodes (both eXeMs and XeMs) con- 
structed with the Planck- ACTPol simulated noise. The 
results from the two sets of eigenmodes are very similar. 
For both XeMs and eXeMs, one mode was sufficient to 
remove the bias {hi = -0.22±0.06 and fii = -0.20±0.06 
respectively). This happens to be in agreement with the 
cutoff mode for the XeMs while with the eXeMs the sec- 
ond mode should also be included. The lower number 
of modes required for the Planck-ACTPol-like case com- 
pared to the ideal experiment is expected due to higher 
sensitivity of the data in the latter to deviations from the 
underlying history. We also tried three modes, with 
no significant detection of the new modes, while the error 
on the XeM 1 increases by a factor of 2. 

4.3. Trajectories 

In this section we investigate the reconstruction of the 
Xc-perturbations using the simulated data to illustrate 
the corresponding uncertainty at different redshifts. The 
left plot in Fig. shows the redshift interval covered by 
500 (51n ATg-trajectories corresponding to an ideal obser- 
vational case with eight XeM included, for the CT2010 
model. The color indicates the number of trajectories 
passing through each (z, 6Xc/Xc) bin, normalized to one 
at each redshift snapshot. The trajectories clearly show 
deviations from the SRS, slowly morphing into the cor- 
rection obtained by CT2010 (the cyan curve). However, 
the recovery is not perfect, as the model of CT2010 has 
non-zero (and relatively significant) projection on higher 
XeMs which are not well constrained by data, and there- 
fore were not included into the analysis. Most obviously, 
corrections to helium recombination are not captured 
well when using only the first few XeMs. These tra- 
jectories do not recover the analytical projection of the 
CT2010 corrections on the first eight XeMs very well 
either. The reason, as discussed before, is that the corre- 
lation of the XeMs induced by the standard parameters 
draws some contribution from the higher absent modes 
which biases the measurement of the first few XeMs in- 
cluded in the measurement. 

To test this impact of higher, excluded modes on the re- 
covered (low XeM) trajectories, we ran simulations with 
the data that only accounted for the contributions from 
these low modes. As expected, in the absence of higher 
modes in the data, the measured XeMs were non-biased 
and thus the highest probability region of the trajectories 
covered the S In X^ curve of the input model. 

Although our basic target is ATc-recovery, the rele- 
vant space for determining how well we have done is 
that of the CMB data, reduced to the power spectra. 



Cj^ and Cf ^. The central and right panels of Fig. [22 
show the dCf^/Cf"^ and 5Cf^/Cf^ trajectories, where 

iCfjCf = (Cf -Cf '*^'*)/Cf '^'^ and Cj'^ is the fiducial 
power spectrum without any perturbations. The trans- 
formation from Xc trajectories to Ci trajectories shows 
a much tighter band around the input signal. This is a 
visual confirmation of the point that some features in the 
(51n ATc which make the trajectories thick do not leave 
a measurable imprint on the 's. Notice there are small 
residual oscillations. They coincide with the peaks and 
troughs of the Ci curves for both TT and EE (which is 
out of phase with TT). One source for the oscillations 
seems to be the eigen mode truncation, as we will see 
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Fig. 22. — Left: The SX^/Xc as measured by a CVL experiment by including eight XeMs (and six standard parameters) in the analysis. 
The colors show the density of trajectories going through each point in the z-5Xc/Xc space, normalized to one at each z. The maximum 
and 1 and 2a widths of the Thomson visibility function have been marked at the bottom of the plot. As this plot and the next ones 
indicate, the main recovery of Xe is the slope of the cure around this visibility peak. Middle and right: similar to the left figure, but for 



oCp /G« trajectories. 
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Fig. 23.— Similar to Fig.p2^but with the first ten eXeMs. 
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Fig. 24. — Similar to Fig. |22| but with the first twenty eXeMs. As this figure demonstrates, including a higher number of modes does not 
necessarily lead to better Ac recovery. Here the recovered Jfc becomes noisier compared to the case with only ten modes included, while 
the On trajectories do not change significantly except for the diminished oscillations around the input model, as discussed in the text. 



later. Using only a limited number of the modes in the 
analysis causes the non-Xe cosmic parameters to try to 
match the injected Xe perturbations. There is also an 
issue of accuracy of the Ci code for some of the distor- 
tions. 

Fig. |23] similarly shows the 2D histograms of trajec- 
tories tor the case with the first ten eXeMs included. 
Around the maximum of the Thomson visibility func- 
tion the Xa reconstruction is slightly stronger and less 
fuzzy than in Fig. [22j with part of the helium recom- 
bination correction being recovered. The improvement 
in the reconstruction is because for the computation of 
these eigenmodes their correlation with the standard pa- 
rameters have been optimally taken care of. In contrast, 
the XeMs used in the previous case are non-optimal if 



no strong additional priors can be placed on the stan- 
dard parameters, leading to confusion i n th e errors and 
the rank-ordering of the modes. Fig. 24 constructed 
with 20 eXeMs included, shows that the oscillation effect 
mentioned above around the input Cg signal is dimin- 
ished (and also partially swamped by the slightly higher 
dispersion around the input curve) for the 20-mode com- 
pared to the ten-mode case. We also see that including 
a higher number of modes does not necessarily lead to 
better recovery. 

Similar trajectories for a Planck-ACTPol-like experi- 
ment are shown in Figs. [25] and [26j with three XeMs 
and three eXeMs as the eigenmodes respectively. The 
trajectories for the XeMs are more widely spread and 
blurred due to experimental noise. The eXeMs perform 
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Fig. 25. — Similar to Fig. |22|but for a Planck- ACTPol-like experiment and with only three XeMs taken into account. 
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Fig. 26. — Similar to Fig. |23|but for a Planck- ACTPol-like experiment and with only three eXeMs taken into account. 



slightly better. However, the overall reconstruction is 
clearly lacking detailed agreement with the full recom- 
bination correction of CT2010. In particular, most of 
the modification during helium recombination is not cap- 
tured, as the corresponding signals can only be picked up 
with higher modes, which in the considered case are not 
constrainable at a significant level. In the 5C^ / Ci plots of 
these Planck- ACTPol-like cases, there is a small disagree- 
ment at high multipoles between the theoretical curve 
and the highest probability region of the chains. That is 
mainly due to mode truncation at a relatively low mode 
number, i.e., three. We tested this by including eight 
modes and as expected, observed a wider spread around 
the input signal with the disagreement diminished. 

Although we do not pl ot the equivalent SCe/Cg for the 
DM case discussed in § |3.5.2[ very similar plots result, 
namely good recovery of the power spectra with a dis- 
persion around the input perturbation signal. 

4.4. Beyond small perturbations 

In this paper it was explicitly assumed that the model 
best explaining the ionization fraction (or the true model 
underlying the ionization history) is only slightly differ- 
ent from our fiducial model, justifying our choice of pa- 
rameter dlnXc- Therefore, the eigenmodes constructed 
for the fiducial model are also very close to the eigen- 
modes for the perturbations to the true X^, the cor- 
rections to the eigenmodes arising from the difference 
between the fiducial and true X^. model being only of 
second order. Under this assumption, a one step search 
for the best-fit parameters suffices to extract the avail- 
able relevant information from the data, provided that 
the minimum required number of modes are included in 
the analysis. Finding the minimum number of required 
modes can by itself involve several parameter estimation 



steps in parameter spaces with different dimensions, the 
criterion being that the best fit values for the standard 
parameters stop changing. That is what was done in 



the examples in this work (§ 4.2), to illustrate how the 
method works. 

However, if the fiducial model is very far from the true 
Xc history, such that the eigenmodes are expected to 
be affected at a significant level, an iterative approach 
toward finding the best modes with their associated am- 
plitudes and errors is required: starting with our best 
guess for the fiducial model, we estimate its deviation 
from the true ionization history using the dataset avail- 
able and the eigenmodes constructed based on this fidu- 
cial model. We then update the model by adding to it 
the measured deviations in the eigenmodes (and the stan- 
dard parameters, if required). This process is repeated 
until the convergence of the model and its eigenmodes. 

However, current constraints seem to indicate that 
such an iterative procedure will not be necessary within 



the standard p icture. For example, as shown by Shaw fc| 
Chluba (2011|, the recombination corrections of CT2011 



are readily incorporated using one calibrated redshift de- 
pendent correction f unction relative to the original re- 
combination model of Seager et al. ( 1999 ). Even for CVL 
errors a second update of the correction functions leads 
to minor effects. Nevertheless, if something more surpris- 
ing occurred during recombination, an iterative approach 
might be required. 

5. CONCLUSION AND DISCUSSION 

CMB data today are becoming so precise that small 
modifications in standard ionization history are impor- 
tant. This impressive progress not only implies that mea- 
surements of the main cosmological parameters are be- 
coming increasingly accurate, but also means that re- 
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maining uncertainties in the recombination dynamics, 
e.g., caused by neglected standard or non-standard phys- 
ical processes, should be quantified. In this work we dis- 
cuss a novel approach to constrain this remaining ambi- 
guity with future CMB data. We performed a principal 
component analysis to find parameter eigenmodes that 
can be used to describe uncertainties in the ionization 
fraction. We constructed Xf. eigenmodes over the red- 
shift range of [200,3000], performing several consistency 
checks to prove the correctness of our method. This ap- 
proach automatically delivers a hierarchy of mode func- 
tions that can be selected according to their error and 
then are added to the standard cosmological parameters 
when analyzing CMB data. 

Due to the strong CMB signal imprinted by hydrogen 
recombination, the most constrained modes are mainly 
localized around z ^ 1100, with some extensio ns t o lower 
and higher redshift regions (see Fig. |4] and [12]). This 
emphasizes that CMB data are very sensitive to small 
changes during hydrogen recombination, while details of 
helium recombination or small changes in the freeze-out 
tail are hard to constrain, unless strong priors on the 
reliability of the hydrogen recombination model are im- 
posed. With the method described here it is possible to 
construct mode functions for different experimental situ- 
ations, also folding in prior knowledge on the recombina- 
tion history using appropriate weight functions and fidu- 
cial Xc models. For example, if there are physically mo- 
tivated and experimentally supported hints toward (sig- 
nificant) changes in the freeze-out tail of recombination, 
e.g., due to energy injection from dark matter annihi- 
lation, we propose a parametrization which weights the 
low redshift part more strongly (see Fig. |9| . 

After we completed th is work, we received a preprint 
(Finkbeiner et al. 2011) which investigated the use of 
CMB data to constrain details of energy injection scenar- 
ios related to decaying or annihilating particles. They 
also used parameter eigenmodes, but these were con- 
structed based on an energy release history which is in 
our language akin to the imposition of a strong prior 
on the recombination dynamics around z ~ 1100 and a 
focus on the freeze-out tail of the recombination. 

We applied the method to different simulated datasets 
with the aim to assess how well future CMB exper- 
iments will be able to constrain modifications to the 
standard recombination scenario. (Current WMAP plus 
ACT and SPT data will provide only relatively weak con- 
straints, but Planck plus ACTPol and SPTPol will con- 
siderably improve the situation.) As a working example 
we used the refined recombination calculations of Cos- 
moRec. For simulated CMB datasets corresponding to 
Planck- ACTPol-like experiments we found that the first 
3 eigenmodes can be rather well constrained. The ad- 
dition of these modes allows us to partially morph from 
the old Recfast X^ model to the new recombination 
history computed with CosmoRec without actually di- 
rectly using the recombination corrections in the analy- 
sis. However, because the first few mode functions are 
strongly localized around z ~ 1100 details during helium 
reco mbin ation and in the freeze-out tail are not captured 
(Fig. 26 ) . The addition of the first 3 eigenmodes is suf- 



up to 10 modes might be constrainable, in this case al- 
lowing us to pick up part of the details during helium 
recombination (Fig. 



23). 



The significance ot the detection of any perturbation 
obviously depends on the underlying ionization history 
of the real data. In the specific CosmoRec example 
for Planck- ACTPol-like experiments, all three modes but 
the first one are consistent with zero. A significant source 
for large errors on the eigenmodes is their correlation 
with the standard parameters. If tight constraints are 
imposed on the standard parameters by non-CMB ex- 
periments such as BAO or SN data, the errors on the 
eigenmodes will be correspondingly reduced. Comparing 
the first rows of Table fl] (where all standard parameters 
are held fixed) and Table Is] (where all standard parame- 
ters are being marginalizea over) illustrates the effect of 
this correlation in the extremes. 

This also shows how important one's knowledge in 
how well elements of recombination are known, expressed 
through prior probabilities, will be. If the uncertainty 
in the ionization history during hydrogen recombination 
can be reliably reduced by other methods then the sensi- 
tivity to small perturbations at higher or lower redshifts 
is enhanced. We note that measurements of the cosmo- 



logical recombination radiation (e.g., see Chluba & Sun- 



yaev 



2006a Sunyaev & Chluba 2009) could in principle 

provide an alternative way ot constraining the recom- 
bination dynamics in the future. In particular, the re- 
combination radiation could exhibit significant features 
if something more unexpected occurred during different 
cosmological epochs (e.g. 
Chluba||2010 ). 



see 



Chluba & Sunyaev 2009 



ficient to remove the biases in the standard parameters, 
however, at the cost of increased error bars. We also 
show that for CVL limited experiments up to / ~ 3500 



In principle, for a complete study of ionization history, 
late reionization should also be included in the analy- 
sis. Ambiguities in the low redshift part of the ioniza- 
tion history may affect the measurements of high redshift 
perturbations and vice versa. However, the main signal 
from the reionization epoch is measurable from the very 
large scale CMB polarization, and the high redshift per- 
turbations of Xc affect anisotropics with smaller angular 
scales. Therefore the signals from these two regions are 
rather uncorrelated. A more complete analysis for the 
whole ionization history or where different parts of it are 
considered simultaneously is for future work. 

An aspect requiring decision when analyzing real data 
is the choice of parametrization. For most of this work 
we weighted the perturbations in X^. by the fiducial his- 
tory. If, for example, the recovered perturbations point 
towards modifications in the freeze-out tail of recombi- 
nation, or if there is strong belief that no sign of sig- 
nificant deviations around the maximum of visibility are 
present, an alternative parametrization which allows bet- 
ter reconstruction of the tail can be constructed, using 
appropriate weight functions that quantify our belief in 
the underlying fiducial model. 

As discussed in § |4.4[ our semi-blind XeMs are de- 
signed to only probe small perturbations about the fidu- 
cial model X^'^. When it comes to real CMB data anal- 
ysis, iterations of X^^ are required to ensure no left- 
over bias remains. We look forward to the application of 
iteratively-improved eigenmodes to the coming high res- 
olution CMB data from Planck, ACTPol and SPTPol. 
We thank Doug Finkbeiner, Tongyan Lin and Olivier 
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